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4.  ABSTRACT  (Continued) 

In  the  case  of  Task  1  we  show  that  supercooled  liquid  water  in 
winter  orographic  clouds  varies  with  the  passage  of  cyclonic  storms,  the 
depth  of  the  cloud  system,  and  proximity  to  mountain  slopes.  Liquid 
water  is  also  found  near  the  tops  of  the  cloud  system  and  at 
temperatures  between  -10  and  -12°C  where  liquid  water  depletion 
processes  are  not  very  efficient.  There  is  also  evidence  of  multiple 
layering  of  liquid  water  ir.  orographic  clouds. 

Sensitivity  studies  with  the  RAMS  cloud  microphysical  module 
reveals  the  extreme  nonlinearity  of  the  microphysical  system.  A 
seemingly  innocent  adjustment  in  crystal  concentration  prediction  can 
lead  to  a  large  response  in  aggregation  rates.  Large  changes  in 
aggregation  rates,  in  turn,  affect  liquid  water  content  since  liquid 
water  is  the  residual  between  production  of  condensate  by  wet  adiabatic 
ascent  (as  radiational  cooling)  and  loss  by  vapor  deposition  and 
accretional  growth  of  ice  crystals.  Aggregation  plays  an  indirect, 
though  important,  role  in  tipping  that  balance.  We  have  mace 
substantial  leaps  toward  the  understanding  of  what  factors  and 
parameterizations  are  of  significance  to  the  prediction  of  supercooled 
liquid  water  content.  At  this  time,  however,  we  do  not  have  the 
capability  of  skillfully  predicting  supercooled  liquid  water  in  a 
variety  of  cloud  situations.  We  also  note  that  it  may  be  necessary  to 
predict  liquid  water  ir.  the  size-range  50-100  in  order  to  predict 
aircraft  icing  conditions  reliably. 

In  the  case  of  Task  2,  a  ID  higher  order  closure  model  has  been 
developed  and  we  have  reviewed  the  literature  regarding  the  properties 
of  upper-level  clouds.  A  major  focus  ir,  the  recent  work  has  been  or. 
formulating  a  new  radiation  parameterization  suitable  for  ice-phase, 
high-level  clouds.  The  microphysical  module  has  also  been  modified  to 
accomodate  the  unique  properties  of  cirrus  clouds.  As  a  first  step  in 
this  direction  the  microphysical  model  has  been  formulated  in  a  unified 
approach  which  allows  the  use  of  a  number  of  different  basic  functions 
for  the  size-distributions  of  drops  and  ice  particles  and  well  as 
predict  the  concentrations  of  hydrometeor  species. 

In  the  case  of  Task  3,  we  have  identified  two  convective  storm 
downdraft  branches  which  contribute  to  strong  low-level  shears  ar.c 
turbulence.  One  branch  called  the  mid-level  branch  is  predominant 
during  the  mature  storm  phase.  The  strength  of  the  latter  branch  is 
controlled  by  the  stability  of  the  atmospheric  boundary  layer  air  ahead 
of  the  storm.  Other  factors,  such  as  precipitation  size  and  phase, 
exert  a  strong  control  on  downdraft  strengths.  In  seme  environments, 
melting  can  account  for  60%  of  the  cooling  of  downdraft  air.  We  have 
also  found  that  as  much  as  30%  of  all  MCC  events  proauce  severe  winds 
with  the  strongest  events  exceeding  45  ms-^  . 

In  Task  4,  we  describe  major  charges  ir.  “-he  CSU  RAMS.  Some  of  the 
most  significant  changes  include  implementation  of  flexible 
interactive-nested  grid  scheme  allowing  an  indefinite  number  of  nests, 
implementation  of  hydrostatic  option,  revis-’er.  of  the  pi  epr„ussui  iuu 
Fortran  77  standard  allowing  greater  model  portability,  and  alteration 
of  the  vectorization  permitting  faster  execution  time  on  computers 
requiring  long  vectors.  The  application  of  the  hydrostatic  version  of 
RAMS  to  mesoscale  convective  complexes  is  also  described. 
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APPENDICES 


1.0  INTRODUCTION 


Ir.  this  report  we  summarize  research  accomplishments  supported  by 
AFC-L  under  contract  #  FI 9628-84-C-0005 .  The  report  summarizes  research 
accomplishments  in  the  following  areas: 

9  Task  1:  A  2D  and/or  3T  model  study  to  define  potential  for  aircraft 
icing . 

»  Task  2:  A  ID  version  of  the  turbulence  closure  model  and  development 
of  a  cirrus  cloud  model. 

ft  Task  3:  Forecasting  and  diagnosis  of  extreme  wind shear  and  turbulence 
associated  with  thunderstorms . 

9  Task  4:  furinary  of  CSU  PAK5  development. 

During  the  course  of  this  contract,  emphasis  was  shifted  from  one 
i  ask  tc  another  ir.  response  to  charges  in  Air  Force  needs,  changes  ir. 
project  staff,  and  charges  ir  staffing  resulting  from  a  slow-down  ir. 
project  spending.  Initially  Task  1  was  perceived  as  a  high,  priority 
:  *  er.  but  as  we  progressed  through,  our  research,  this  area  was  perceived 
tc  te  of  lesser  importance  to  the  Air  Force  than  ether  tasks.  Likewise 
Task  3  was  a  major  focus  cf  this  research,  during  Kevin  Kr.upp’s 
disrertat jor.  rtudy.  5ir.ee  completing  his  thesis  this  task  has  received 
lesser  emphasis,  though,  our  analysis  of  Perecho-type  events  associated 
with  h'CCs  suggests  that  it  should  remain  ar.  area  of  concern  to  the  Air 
Force.  Task  2  has  undergone  considerable  evolution  during  this  project. 
It  initially  started  out  with  the  intent  of  developing  a  simplified  ID 
model  to  predict  cumulus  cloud  cover  and  middle  and  high  clouds.  The 
erphasis  has  shifted  miore  to  the  physics  (radiative  properties  and 
micrcphysics )  of  cirrus-type  clouds  and  middle-level  stratiform  cloud 
using  PAIN?  the  dynamic  model.  The  shift  in  emphasis  in  Task  2  was  ir. 


response  to  a  perceived  change  1  r.  Air  Force  interests  and  also  the 
change  in  focus  of  Piotr  Flatau's  research  interests.  We  hope  tc 
continue  a  focus  in  Task  2  under  Air  Force  support.  Task  4  is  only 
briefly  summarized  here  as  the  objective  of  this  task  was  to  keep  AFGL 
personnel  informed  about  changes  and  developments  in  RAMS  on  a  regular 
basis.  This  has  been  accomplished  through  regular-  quarterly  reports  arc 
personal  discussions. 

2.0  MODEL  AND  OBSERVATIONAL  STUDIES  RELATFD  TO  DEFINING  POTENTIAL 

AIRCRAFT  ICING  CONDITIONS 
2.1  Introduction 

As  part  cf  several  ongoing  research  projects,  using  the  RAMS 
modelling  system  as  ar.  investigation  tool,  we  have  studied  the 
microphysical  development  cf  several  different  classes  cf  cloud  forms. 

As  a  result,  we  have  beer,  led  tc  greatly  improve  cur  bulk  microphysics 
parameterization  and  become  aware  cf  some  of  the  more  fundamental  issues 
concerning  aircraft  icing.  In  this  section,  we  review  the  progress  we 
have  mace . 

Our  test  bed  for  cicrophysic s  studies  are  observed  orographic 
clouds.  In  particular,  we  have  used  results  of  the  COSE  (Colorado 
Orographic  Seeding  Experiment)  rur.  by  Dr.  L.  C.  Grant,  of  Colorado  State 
University.  This  ongoing  observational  program  has  provided  us  with  a 
substantial  data  base  against  which  we  can  test  our  r.icr ophysics 
hypotheses.  The  overwhelming  advantage  cf  studying  orographic  clouds  :: 
that  their  dynamics  are  relatively  stationary  and  deterministic;  a 
simple  function  of  the  mean  ambient  wind.  Parcel  history  can  be 
observed  spatially,  varying  as  distance  from  the  orographic  barrier. 

With  minimum  uncertainty,  we  are  thus  able  to  compare  our  model 


pr  edit  tic  nr.  directly  to  observations.  This  task  is  rather  elusive  wher. 
cr.e  considers  the  case  of  non-steady  convective  clouds. 

2.2  Review  of  Previous  Research 

The  liquic  water  content  at  any  location  ir.  a  clcud  is  the 
integrated  consequence  of  several  liquid  water  production  prccesses  l FP ) 
and  liquid  water  depletion  processes  (FD).  Liquid  water  production  is 
determined  ty  the  cooling  rate  of  the  air  as  it  rises.  The  cooling 
rate,  in  turn,  is  primarily  due  to  adiabatic  ascent.  Net  radiative 
cooling  is  also  a  factor,  at  least  in  weaker  wind  situations  ir,  the 
mountains  and  ir.  the  more  stratiform  cloud  systems.  The  major  liquid 
water  depletion  prccesses  belcw  OC  are  vapor  depcsitional  growth,  of  ice 
crystals  ard  riming  ty  clcud  droplets.  The  depletion  rates  are 
therefore  a  function  cf  the  concentration  of  ice  particles,  the 
temperature  end  h =  of  the  ice  particles,  the  cloud  supersaturaticr, 
and  liquic  water  content. 

Wner.ever  FF  >  FP,  liquid  water  will  accumulate  ir.  the  cloud  to 
levels  that  car.  te  collected  or.  airframe  surfaces.  If,  however,  FP  < 

FP .  liquid  water  will  be  depleted  so  rapidly  that  only  solid 
precipitation  elements  will  be  detectable.  As  a  simple  example  of  the 
competition  between  FP  and  FD,  consider  the  diagram  shown  ir.  Fig.  2.1 
obtained  by  Chappell  (1970)  for  orographic  clouds.  Tice  figure 
illustrates  estimated  water  production  rates  for  airflow  over  the 
Climax,  Colorado  mountain  barrier  as  a  function  of  the  500  mb 
temperature  ir.  that  region.  Chappell  used  the  500  mb  temperature  as  a 
crude  index  of  the  cloud  top  temperature  in  that  region.  As  the  cloud 
top  temperatures  and  also  cloud  bases  (because  the  bases  are  limited  by 
the  terrain  height)  become  colder,  saturation  mixing  ratios  diminish. 
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Fig.  2.1. 


Distribution  of  non-seeded  precipitation  at  HAD  as  a  function 
of  500  r.b  temperature  compared  to  a  theoretical  distribution 
computed  using  the  mean  diffusicnal  model.  Precipitation 
data  are  from  Climax  I  sample  (251)  and  values  are  a  running 
mean  over  a  two-degree  temperature  interval.  [From  Chappell, 
1970.  ] 


thereby  decree s ir.g  FP.  The  lcng-dashed  curve  illustrates  Chappell's 
estimate  of  FD.  Chappell  assur.ed  that  the  major  depletion  process  was 
vapor  deposition  growth  of  ice  particles.  He  also  assumed  that  the 
oor.cer.traticr.  of  ice  crystals  could  be  determined  by  the  Fletcher  ice 
nucleus  spectrum.  This  figure  illustrates  that  at  temperatures  warmer 
than  -20  c,  FP  >>  FD,  thus  the  opportunities  for  significant  liquid 
water  production  are  large.  At  temperatures  colder  than  -20  C,  the  more 
numerous  ice  crystals  are  able  to  deplete  liquid  water  at  a  rate  faster 
thar.  it  is  being  produced.  Also  shown  ir.  Fig.  2.1  are  the  observed 
surface  precipitation  rates.  At  temperatures  colder  than  -20  C,  the 
observed  precipitation  is  comparable  ir.  magnitude  to  the  calculated  FT. 
This  sirply  shews  that  at  colder  temperatures,  FD  is  limited  by  FP  ant 
little  excess  liquid  water  is  produced.  At  temperatures  warmer  than  -20 
C,  the  otsc-ved  precipitation  rate  is  considerably  loss  than  Chappell’s 
estimate  of  FP ,  shewing  an  opportunity  for  production  of  considerable 
liquid  water.  Note  that  the  observed  precipitation  rate  exceeds 
Chappell's  estimate  of  FD.  This  is  probably  because  Chappell  assumed 
only  vapor  deposition  contributed  tc  FD  whereas  mining  growth  of  ice 
crystals  could  have  been  significant  at  these  warmer  temperatures  where 
high  liquid  water  contents  are  likely.  Alsu  Chappell  used  the  Fletcher 
ice  nucleus  curve  which,  as  discussed  below,  is  frequently  a  poor 
estimate  of  observed  ice  crystal  concentrations  at  warmer  temperatures 
because  it  dees  not  consider  the  effects  of  ice  multiplication. 

The  observed  structure  of  the  supercooled  liquid  water  region  ir 
crographic  clouds  is  quite  variable  ar.d  depends  on  the  particular 
dynamical  properties  of  the  cloud  system.  As  we  shall  see,  the 
detection  of  supercooled  water  also  depends  on  the  particular  water 


T 


V 


sensing  systetr. . 

Stable  orograhic  clouds 

The  simplest  cloud  system  ir.  which  to  interpret  the  behavior  cf 
supercooled  water,  is  the  "pure"  orographic  cloud  system.  This  clcuc 
system  occurs  over  mountainous  regions  during  periods  when  large-scale 
frontal  activity  is  absent.  They  typically  occur  wren  a  strong  cross¬ 
barrier  pressure  gradient  is  present  driving  strong  barrier- normal  winds 
containing  moist  air.  These  are  the  cloud  systems  that  are  often  well 
simulated  by  linear  orographic  flow  models.  As  found  by  Hobbs  (1973) 
over  the  Cascade  Mountains  of  Washington,  and  Fauter  (1985)  over  the 
Park  Range  of  northern  Colorado,  supercooled  liquid  water  is  largest  ir. 
regions  of  strong  orographic  forcing  or.  the  windward  side  cf  steep  rises 
in  topography.  Raucer  also  noted  that  in  some  cases  liquid  water  was 
observed  to  vary  uniformly  through  the  depth  cf  the  orographic  clcuc 
layer  without  any  discerr.atle  changer  ir.  the  dynamical  properties  cf  the 
cloud  system.  The  only  detectable  change  in  clcuc  structure  was  a 
reduction  ir,  precipitation  rate  which  indicated  that  variations  ir.  ice- 
crystal  concentrations  may  be  responsible.  The  causes  of  such 
concentration  variations  are  r.ct,  well  known,  but  the  evidence  presentee 
by  PeMott  et  al .  (  1986  )  suggests  that  enhanced  r.ucleation  of  :ce; 
particles  occurs  ir.  localized  region;  of  high  super-saturations 
associated  with  upward  motions  rear  cloud  top. 

Supercooled  liquid  water  in  orcgrahicall y-modif ied  cyclonic  storm; 

Variations  in  the  distribution  of  liquid  water  ir.  orographic  clcuc: 
are  frequently  associated  with  variations  in  thermocyr.amic 
stratification,  flow  and  vertical  motion  as  cyclonic  storm  systems  pass 


over  the  mountain  barrier.  Eased  cr  aircraft  observations  cf  liquid 


water  and  rri c-rophy sical  structure  of  orographic  clouds  over  the  Sar.  Juan 
mountains  of  southern  Colorado,  Marwitz  (1980)  and  Cooper  and  Marwitz 
(1980)  identified  four  stages  in  thermodynamic  stratification  of  a 
passing  storm  system.  Unfortunately  they  did  net  identify  there  regimes 
with  frontal  boundaries  or  features  of  the  Norwegian  Cyclone  model. 

This  may  he  due  to  the  fact  that  extra-tropical  cyclones  are  sc  strongly 
altered  in  structure  by  the  terrain  features  that  classical  warm  fronts 
arc'  cole  fronts  could  not  be  unambiguously  identified. 

Nonetheless,  they  identified  four  stages  in  the  evolution  of  the 
thermodynamic  structure  of  the  storm  system  shewn  in  Fig.  2.2,  which  can 
roughly  be  identified  with  the  passage  of  a  cyclonic  stern:.  Fig.  2.3 
illustrates  a  conceptual  model  of  the  liquid  water  distribution  for 
three  cf  the  four  stages.  Turing  the  "stable  stage"  which  generally 
cor  responds  to  the  pre-frontal  stage  cf  a  cyclonic  storm,  moisture 
supplied  by  the  rising  warm  "conveyer  Pelt"  contributes  to  a  deep 
stable  cloud  layer  with  cold  top:  .  The  eicud  system  is  nearly 
completely  glaciated,  owing  to  the  high  ice  crystal  concentrations 
associated  with  the  very  cold  tic  ud  tops.  Liquid  water  was  observed 
only  in  short-lived  patches  where  gravity  waves  penetrate  transient 
regions  of  rapid  ascert.  Any  liquid  water  that  is  produced  is  small, 
being  less  than  0.1  g/m  .  As  can  be  seen  from  Fig.  2.3a,  substantial 
blocking  of  the  lcw-level  flow  is  indicated  during  this  period.  As  a 
cyclonic  storm  advances  across  the  mountain  barrier,  the  atmosphere 
transforms  to  near  neutral  stability  through  much  of  the  depth  of  the 
troposphere .  Liquid  water  contents  rise  to  values  on  the  order  of  0.3 
ftr?  (see  Fig.  2.3b).  Following  the  sassage  of  the  500  mb  level  trough 


and  in  seme  cases  a  discernafclf  front,  the  atmosphere  transformed  to  be 


2.2.  Schematic  height-time  cross  section  for  a  typical  storm  ever 
the  San  Juans.  The  solid  lines  are  lines  of  constant 
equivalent  potential  temperature  (0e) ,  and  the  shaded  area 
indicates  convective  instability.  The  typical  sequence  of 
500  mb  temperature  is  shown  at  the  top.  The  reference 
baseline  is  -23°C.  The  level  below  which  the  wind  direction 
are  <  160c  is  indicated  by  circles  and  is  defined  as  blocked 
flew.  The  four  stages  of  the  storm  and  a  typical  seeding 
period  are  indicated.  [From  Cooper  and  Karwitz,  1980.] 
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Fig.  2.3.  Typical  appearance  of  San  Juan  storms  during  the  (a)  stable, 
(t)  neutral,  and  (c)  unstable  storm  stages.  The  airflow  is 
from  the  left  at  ~  15  m  s~3 .  The  shaded  regions  are  regions 
where  liquid  water  is  generally  found,  and  some  typical 
liquid  water  contents  (g  nT3)  are  indicated.  Solid  lines 
with  arrows  show  air  trajectories,  and  the  irregular  lines 
indicate  cloud  boundaries.  Some  typical  ice  crystal 
trajectories  are  also  shown  as  dashed  lines.  The  farthest 
downwind  trajectory  shows  the  fall  of  a  hydrometeor  with  a 
terminal  velocity  of  1  it  s~3 .  [Frcm  Cooper  and  Marwitz , 
1980.1 


cor.vectively  unstable  at  lower  levels  upwind  of  the  mountains  where  a 


surface  zone  of  convergence  was  observed.  Cooper  ar.d  Marwitz  reported 
that  liquid  water  contents  were  highest  during  this  stage  of  the  storr. 
system,  approaching  1.0  g/m3  in  the  upper  cloud  levels.  Hobbs  (1975) 
and  Heggli  et  al .  (  1983  )  also  found  the  highest  liquid  water  contents  ir. 
post-frontal  convective  clouds  observed  over  the  Cascade  mountains  cf 
the  state  of  Washington  and  the  California  Sierras,  respectively.  Likc- 
Cooper  and  Marwitz,  their  observations  were  primarily  obtained  by 
aircraft  penetrations.  One  deficiency  of  aircraft  observations  is  that 
for  safety  reasons,  aircraft  are  limited  to  flying  at  heights  600  m  tc 
1000  m  above  mountain  tops.  Also  in  the  case  of  convective  cells 
imbedded  in  stratiform  cloud,  time  continuity  of  the  convective  cells 
cannot  be  obtained.  Pecently  remote  probing  systems  such,  as  microwave 
radiometers  (Hogg,  et  al.  .  1983  )  and  polarization  diversity  lidar 
(Cassen,  1984)  are  providing  an  opportunity  to  sense  the  presence  cf 
liquid  water  in  regions  previously  inaccessible  by  aircraft.  We  shall 
see  that  data  obtained  with  these  new  observational  systems  is 
substantially  altering  our  view  of  the  distribution  of  liquid  water  in 
orographic  clouds. 

Using  a  combination  of  aircraft  observations  and  passive  microwave 
radar  observations,  Rauber  (1985),  Rauber  et  al .  (1986)  and  Rauber  anc 
Grant  (1986)  have  described  the  distribution  of  supercooled  liquid  water 
observed  in  the  Park  Range  of  northern  Colorado  during  passage  cf 
extra-tropical  cyclones.  They  noted  that  in  the  region  just  aheaa  cf 
and  during  frontal  passage,  clouds  are  characterized  by  wide  area 
stratiform  regions  with  frequently  superimposed  deeper  clouds  and 
heavier  precipitation.  Many  of  the  heavier  precipitating  regions  are 


convective  but  sene  occur  ir.  a  stably  stratified  or  neutral  atmosphere. 
Liquid  water  production  in  the  more  stratiform,  region  of  the  pre-frontal 
environment  is  strongly  influenced  by  the  orientation  of  the  approaching 
front  relative  to  the  orientation  of  the  mountain  barrier.  In  the  case 
cf  the  Park  Far.ge,  mid-level  flow  with  a  strong  westenly  component 
generates  strong  orographic  lifting  and  is  generally  of  higher  moisture 
content.  Convection  was  often  observed  in  the  pre-frontal  environment. 
It  was  sometimes  organized  into  a  banded  structure.  Liquid  water  was 
found  in  the  convective  clouds  during  their  initial  stages  of 
development.  However,  the  liquid  water  was  rapidly  depleted  by  very 
efficient  precipitation  processes.  In  general,  liquid  water  amounts 
reached  a  minimum  within  15-20  minutes  after  the  onset  of  convection , 

The  post-frontal  clouds  over  the  Park  Range  were  also  characterized 
by  wide  area  stratiform  regions.  Orographic  forcing  is  particularxy 
strong  because  the  typically  strong  westerly  to  northwesterly  mid-level 
winds  are  normal  to  the  barrier  orientation.  The  highest  liquid  water 
amounts  were  found  10-15  km  upwind  of  the  barrier  crest  and  reach  a 
maximum  over  the  windward  slopes  where  the  condensation  supply  rates  are 
at  maximum  due  to  orographic  lifting.  Two  regions  of  convective 
instability  were  observed  in  the  post-frontal  region.  The  first 
occurred  just  after  frontal  passage,  and  the  second  during  the  passage 
of  the  500  mb  level  trough  as  the  system  dissipated.  The  distribution 
of  liquid  water  in  the  post-frontal  zone  was  complicated  by  the  fact 
that  superimposed  on  steady  orographic  production  of  liquid  water  was 
the  mere  transient  convective  activity.  Again  the  convection  produced 
liquid  water  early  in  its  lifecycle  but  rapidly  produced  substantial 
numbers  of  ice  crystals  that  depleted  liquid  water  ir.  the  convective 


regions  as  well  as  the  lcwer-levei  feeder  regions. 

Fig.  2.4  is  a  schematic  illustration  of  the  observed  distribution 
of  supercooled  liquid  water  in  (a)  shallow  orographic  clouds,  (tl  deep 
stratiform  cloud  systems,  and  (c)  convective  regicr,  imbedded  ir.  a 
shallow  stratiform  cloud  system.  In  all  cloud  systems  liquid  water  was 
observed  in  the  "feeder"  zone  directly  upwind  and  over  the  barrier 
crest.  A  second  region  of  liquid  water  was  inferred  between  cloud  base 
and  the  -10  C  level.  In  this  region,  planar  crystals  falling  into  it  do 
not  grow  rapidly  by  vapor  depositicr,  growth.  Also  the  more  rapidly 
growing  needles  and  columnar  form  of  crystals  were  characteristically 
absent.  Thus  the  liquid  water  depletion  mechanisms  (FD)  were  at  a 
minimum,  yet  water  production  (FP)  at  these  warmer  temperatures  wan 
quite  high.  A  third  region  near  the  tops  of  the  clouds  illustrated  ir. 
Fig.  2.4  is  also  a  region  of  frequent  liquid  water  occurrence.  Rauter 
and  Grant  (1986)  argued  that  this  zone  of  liquid  water  is  a  result  of  an 
imbalance  between  liquid  water  production  by  adiabatic  cooling  and 
depletion  by  vapor  deposition.  They  argued  that  large  ice  crystals 
capable  of  extracting  large  quantities  of  vapor  as  they  grow  (and  hence 
depleting  liquid  water  by  the  Eergeror.-Feir.deisen  process)  rapidly 
settle  out  of  the  cloud  top  layer.  Left  behind  are  numerous  small  ice 
crystals  which  are  incapable  of  extracting  water  substance  at  the  rates 
being  supplied  by  adiabatic  ascent.  Another  possible  contribution  to 
the  observed  liquid  water  near  cloud  top  is  that,  liquic  water  production 
could  be  intensified  by  radiative  cooling.  Radiative  cooling  rates 
estimated  with  the  Cher,  and  Cotton  (19  86)  stratocumulus  model  suggest 
that  it  can  approach  200°C/day  in  a  layer  of  5  mb  thickness  near  cloud 
top.  For  moist  adiabatic  ascent,  the  cooling  rate  (CR)  is 
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Fig.  2.4.  Conceptual  models  of  supercooled  liquid  water  distribution  in 
Park  Range  cloud  systems  for  (A)  a  shallow  stratiform  cloud 
system  with  a  warm  cloud  top  temperature  (CTT  >  -22°C);  (B)  a 
deep  stratiform  system  with  a  cold  cloud  top  temperature;  and 
(c)  a  deep  convective  region  embedded  in  a  shallow  stratiform 
system.  The  characteristic  precipitation  type  is  listed  for 
each  case  along  the  bottom  of  the  figure.  Magnitudes  of 
liquid  water  contents  measured  with  aircraft  or  at  SPL  are 
shown  on  the  figure.  fFrom  Pauber  and  Grant,  1986.] 


CR  -  665w  °C/day , 


where  w  is  the  vertical  motion  ir.  meters  per  second.  Thus  ar.  updraft  of 
0.3  m/s  is  required  to  produce  the  same  amount  of  cooling  anc  liquid 
water  as  cloud  top  radiation.  Cloud  top  radiative  coding  car.  therefore 
contribute  to  an  imbalance  between  FP  and  FD.  Radiative  cooling  is 
probably  relatively  more  important  ir.  the  more  stratiform  cloud  region 
that  blankets  the  upwind  side  of  the  barrier. 

Observations  of  supercooled  liquid  water  over  the  Sierra  Nevada 

Aircraft  observations  of  supercooled  liquid  water  in  winter  cloud 
systems  over  the  Sierra  Nevada  (Heggli  et  al . ,  1983;  Lamb  et  al . ,  1976) 
suggested  that  the  largest  amounts  of  supercooled  liquid  water  resided 
ir.  convective  towers.  Fere  recent  c bservatior.s  using  passive  microwave 
radiometers  as  well  as  observations  cf  heavily  rimed  ice  crystals  from 
shallow  orographic  clouds  led  Reynolds  and  Dennis  (1986)  to  concluae 
that  the  shallow  nor.-cor.vective  clcucs  contain  substantial  amounts  cf 
supercooled  liquid  water.  Heggli  (  1986)  noted  that  about  75®:  of  the 
time  saturated  clouds,  possibly  containing  liquid  water,  existed  below 
minimum  aircraft  flight  altitudes.  Heggli  noted  that  the  vertical 
distribution  of  liquid  water  was  fci-nodal.  A  low-level  maximum  was 
observed  between  -2  and  -4  C  and  a  secondary  maximum  similar  to  that 
observed  by  Rauber  and  Grant  (1986)  occurred  between  -K1  and  -12  C. 
Heggli  also  noted  that  supercooled  liquid  water  cic  net  always  c  .cur  ir. 
single  layers.  They  often  observed  that  supercooled  liquid  water 
occurred  in  multiple  layers  which  were  separated  by  100m  or  mere  of 
unsaturated  air.  The  mechanisms  responsible  for  such  layering  of  the 
liquid  water  zones,  however,  have  not  been  identified. 

Heggli  and  Reynolds  (1985)  and  Feynolds  and  Dennis  (1986) 


V 


interpreted  the  behavior  of  liquid  water  variations  during  the  passage 
of  cyclonic  storms  over  the  Sierra  Nevada  in  terms  of  Browning  and 
Monk's  (1982)  split-front  model  shown  in  Fig.  2.5.  Following  the 
passage  of  the  frontal  rainbanc  where  the  heaviest  precipitation  was 
observed,  a  shallow  orographic  cloud,  often  of  the  form  of  a 
stratocumulus  cloud,  prevails.  It  is  ir,  the  shallow  orographic  cloud 
layer  that  the  highest  sustained  liquid  water  contents  are  observed. 
Reynolds  (1986)  attributed  this  rise  in  supercooled  liquid  water  to  the 
disappearance  of  ice  crystals  settling  from  higher  clouds.  In  terms  of 
the  " seeder- feeder "  mechanism,  the  deep  clcuds  which  prevail  during 
the  passage  of  the  elevated  frontal  band  yield  an  environment  that  is 
very  favorable  for  the  operation  of  the  seeder-feeder  process.  High- 
level  feeder  crystals  readily  deplete  liquid  water  being  generated  by 
the  lower-level,  imbedded  orographic  feeder  clouds.  Following  the 
passage  of  the  frontal  rainbar.d,  only  shallow  cloud  layers  prevail  thus 
essentially  shutting  down  the  seeder-feeder  process.  This  is  consistent 
with  Kuciauskas’  (1986)  observation  that  substantial  supercooled  liquid 
water  content  amounts  occurred  several  hours  following  the  passage  cf 
the  elevated  cloud  layer. 

In  summary,  we  have  seen  ir.  this  section,  that  supercooled  liquid 
water  in  winter  orographic  clouds  is  quite  variable  in  its  spatial  and 
temporal  distribution.  It  varies  with  the  passage  of  cyclonic  storm 
systems,  the  depth  of  the  cloud  system,  and  proximity  to  mountain 
slopes.  Liquid  water  is  also  found  near  the  tops  of  the  cloud  system 
and  at  temperatures  between  -10  and  -12  C  where  liquid  water  depletion 
processes  are  not  very  efficient.  There  is  also  evidence  of  multiple 
layering  of  liquid  water  in  some  orographic  clouds.  The  mechanisms 
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Fig.  2.5.  Schematic  portrayal  cf  a  split-front  with  the  warm  conveyer 

belt  undergoing  forward  sloping  ascent,  but  drawing  attention 
to  the  split  front  characteristic  and  the  overall 
precipitation  distribution:  (a)  plan  view,  (b)  vertical 
section  along  AB  in  (a).  In  (a)  t'U  represents  the  upper  cole 
front.  The  hatched  shading  along  UU  and  ahead  of  the  warm 
front  represents  precipitation  associated  with  the  upper  cold 
front  and  warm  front,  respectively.  Numbers  in  (t)  represent 
precipitation  type  as  follows:  1,  warm  frontal 
precipitation;  2,  convective  precipitation-generating  cells 
associated  with  the  upper  cold  front;  3,  precipitation  for 
the  upper  cold  frontal  convection  descending  through  an  area 
of  warm  advection;  4,  shallow  moist  zone  between  the  upper 
and  surface  cold  fronts  characterized  by  warm  advection  and 
scattered  outbreaks  of  mainly  light  rain  and  drizzle;  5, 
shallow  precipitation  at  the  surface  cold  front  itself. 
[Adapted  from  Browning  and  Monk,  1982,  Browning,  1985,  and 
Reynolds  and  Dennis,  1986.] 


responsible  for  such  layering  of  liquid  water  are  not  known  at  the 
present  time,  however. 

2.3  RAMS  Model  Development  Appropriate  to  Aircraft  Icing  Studies 

Ir.  the  past  several  years,  we  have  centered  our  mierophysical  work 
around  the  investigation  of  the  aggregation  process.  Rauber  (1985) 
concluded  in  his  study  of  COSF  data  that  aggregation  was  a  highly 
misunderstood  but  extremely  important  phenomena  which  greatly  modulates 
the  precipitation  process.  Our  work  investigating  MCCs  (Mesoscale 
Convective  Complexes)  also  reveals  a  great  importance  of  aggregation  in 
extensive  stratiform  anvil-like  clouds  to  the  overall  precipitation 
process  (Yeh  et  al.,  1987). 

Our  need  to  model  the  aggregation  process  led  us  to  develop  a 
crystal  concentration  and  aggregation  prediction  scheme  published  by 
Cotton  et  al .  (1986).  Our  attempt  to  adjust  the  prediction  scheme  tc 
match  observations  proved  difficult.  As  part  of  this  investigation,  we 
participated  in  the  International  Cloud  Modeling  Workshop,  held  in 
Irsee,  F.G.F.  ir.  1985.  As  part  of  that  workshop  we  attempted  to 
simulate  a  "blind"  observed  orographic  seeding  case  in  order  to 
demonstrate  the  true  skill  of  our  prediction  scheme.  Rauber  selected  a 
blind  case  where  aggregates  did  not  form  in  order  to  see  if  we  could 
properly  sense  the  conditions  for  aggregation.  Although  we  were  able  to 
adjust  the  model  to  match  an  open  case  study,  our  attempts  at  the  Blind 
case  study  met  with  failure  ir,  that  aggregates  were  predicted.  At  this 
time,  we  are  attempting  to  learn  more  about  aggregation  so  that  we  can 
build  the  needed  skill  into  our  parameterization. 

This  ongoing  study  has  illustrated  the  great  non-linearity  of  this 
microphysical  process.  A  seemingly  innocent  adjustment  in  crystal 


concentration  prediction  can  lead  to  a  large  response  in  aggregation. 

The  reactions  to  the  changes  we  make  are  seldom  intuitively  obvious 
until  investigated  in  detail.  As  was  shown  in  section  2.2,  liquid  water 
content  is  a  residual  between  production  of  condensate  cooling  of  moist 
air  and  loss  by  depositional  and  accretional  growth  of  ice  crystals. 
Aggregation  plays  an  indirect  but  significant  role  in  tipping  that 
balance.  We  find  that  aggregation  rapidly  depletes  the  number  of  ice 
crystals  so  that  despite  increased  riming,  less  depositional  ice  growth 
leads  to  the  existence  of  greater  liquid  water.  From  observations, 
Rauber  (1985)  demonstrated  that  liquid  water  tends  to  form  at  cloud  top. 
He  hypothesized  this  may  be  due  to  the  depletion  of  ice  crystals  by 
precipitation  at  cloud  top  again  tipping  the  balance  toward  growth  of 
cloud  water.  More  recently.  Cotter,  (see  Cotton  and  Anth.es,  1987) 
suggested  that  longwave  radiation  divergence  at  cloud  top  may  also 
assist  in  tipping  the  balance  toward  production  of  liquid  water  there. 
This  is  currently  being  investigated. 

Our  recent  observational  studies  of  MCC  stratiform  rain  have 
suggested  that  aggregation  may  be  occurring  far  core  frequently  than  we 
have  yet  realized.  Yeh  et  al.  (1987)  have  found  a  significant 
aggregation  process  at  temperatures  colder  than  -12C.  Passarelli  (1978) 
and  later  Rauber  (1985)  suggested  mechanical  interaction  of  dendrite 
crystals  to  be  the  dominant  production  mechanism,  confining  the  process 
to  the  dentrltic  zone  between  -12  and  -15  C.  However,  Rauber  pointed 
out  an  important  contribution  by  spatial  dendrites  at  temperatures  below 
-12  C  in  initiating  the  process. 


We  have  become  convinced  that  aggregation  is  a  key  element  needed 
to  predict  all  facets  of  microphysical  structure.  Ferhaps  supercooled 
liquid  water  content  is  the  most  sensitive  to  its  prediction.  The 
prediction  of  aggregation  itself  is  highly  sensitive  to  two  primary 
factors.  The  first  is  the  aggregation  efficiency.  The  second  is  the 
distribution  and  number  of  crystals.  This  brings  the  ice  crystal 
nucleation  into  the  problem  in  a  big  way.  Our  best  guess  on 
sorption/depositional  nucleation  rates  is  an  assumed  movement  along  the 
empirically-derived  Fletcher  curve,  which  is  only  a  function  of 
temperature.  We  know  little  about  how  the  nucleation  rates  vary  with 
pressure,  altitude,  air  mass,  supersaturation,  etc.  Grant  (personal 
communication)  has  pointed  out  that  ice  nucleus  measurements  during  the 
CLIMAX  experiment  generally  departed  significantly  from  the  Fletcher 
curve,  being  substantially  less  at  temperatures  cclder  than  -20  C. 
Contact  nucleation  at  least  using  Young’s  (1974)  model  also  does  not 
appear  to  be  adequate.  Thus,  ere  of  the  central  most  important 
parameters  to  the  aggregation  process,  i.e.  crystal  production  by 
nucleation,  still  cannot  be  predicted  within  several  orders  of  magnitude 
routinely.  We  suggest  that  more  resources  be  placed  toward  the 
formation  of  theoretical  and  laboratory  studies  of  models  for  all 
microphysical  processes  and,  ir,  particular,  primary  and  secondary  ice 
crystal  production. 

In  studying  our  crystal  concentration  parameterization,  we  have 
selected  several  major  problems  in  need  of  immediate  attention.  The 
most  important  is  that  we  assume  a  single  category  for  all  pristine 
crystals.  Because  of  the  non-lir.ear  production  rates  with  temperature 
suggested  by  the  Fletcher  curve,  a  numerical  model  rapidly  becomes 
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dominated  by  the  large  numbers  of  small  ice  crystals  initiated  at 
temperatures  colder  than  -25  C.  In  an  updraft  this  essentially  shrinks 
crystals  which  should  be  growing,  ultimately  inhibiting  precipitaticr. 
sedimentation  back  into  lower  levels.  We  have  decided  to  improve  the 
parameterization  by  further  subdividing  pristine  crystals  into  snow-like 
crystals  of  relatively  large  size  and  terminal  velocity  and  very  small 
pristine  crystals.  This  allows  the  coexistence  of  small  and  large 
crystals  at  cold  temperatures,  enabling  the  continuous  growth  of 
crystals  in  the  updraft.  Also,  instead  of  a  constant  size  distribution, 
we  will  assume  a  variable  size  distribution  for  the  snow.  Hence,  we  are 
new  planning  to  consider  crystals  in  their  pristine  form,  mature- 
precipitating  form,  aggregated  forms  as  well  as  graupei  particles. 

At  this  pcir.t  in  the  development  of  our  microphysical  model,  we 
have  found  that  the  formulation  has  matured  to  the  point  where  we  are  in 
need  of  a  totally  new  approach  tc  the  code  emphasizing  user- 
friendliness,  flexibility,  and  parallels  between  all  of  the 
formulations.  We  also  found  a  need  for  transportability  so  that  the 
code  can  be  used  with  other  dynamical  models  outside  of  our  system,  '..'t 
are  therefore  now  completing  the  development  cf  a  new  generalized 
microphy sical  code  capable  of  quick  and  significant  changes  ir,  the 
scheme,  and  incorporating  the  newer  physics  described  above.  This 
parameterization  will  be  defaulted  to  the  original  parameterization 
reported  by  Cotton  et  al .  (1986),  but  is  capable  of  instant 
modifications  by  the  addition  of  new  hydrometeor  categories,  alternate 
size  distribution  assumptions,  alternate  diagnostic  and  prognostic 
sememes  for  distribution  parameters,  and  other  physics  modules.  This 


code  became  operational  during  the  fall  of  1987.  We  expect  that  it  will 
greatly  accelerate  our  rate  of  development. 

Ir.  summary,  our  work  in  microphysics  over  the  past  several  years 
has  demonstrated  that  the  existence  of  supercooled  water  and  thus  the 
potential  for  aircraft  icing  is  perhaps  one  of  the  most  difficult  of  all 
micrcphysical  paramters  to  predict,  primarily  because  it  tends  to  exist 
as  a  residual  of  large  rates  of  water  production  and  large  rates  of 
precip itaticr.  removal  processes.  We  have  made  substantial  leaps  toward 
the  understanding  of  what  factors  are  of  significance .  At  this  time,  we 
do  not  have  the  capability  to  skillfully  predict  supercooled  water, 
however  we  are  continuing  to  work,  both  tnrough  observation  and 
modelling,  to  achieve  that  goal.  We  also  recognize  that  prediction  of 
the  amount  of  supercooled  water  is  not  the  only  requirement  for  aircraft 
icing  conditions.  The  droplet  spectra,  particularly  the  presence  of 
droplets  ir.  the  size  range  50-100  pm  greatly  aggrevates  aircraft  icing 
conditions.  'i’r.for tur.ately,  this  is  precisely  the  size  range  that  our 
rccei  has  the  least  skill  ir.  predicting  the  amounts  of  supercooled 
water.  Another  approach  such  as  the  multiple  log-normal  basis  function 
approach  of  Clark  and  Hal  1  (  1  983  )  may  be  required. 

3.0  A  ID  VERSION  OF  THE  TURBULENCE  CLOSURE  MODEL  AND  DEVELOPMENT  OF  A 

CIRRUS  CLOUD  MODEL 
3  .1  Introduction 

Research  in  this  section  began  as  an  extension  of  the  Chen  and 
Cottor.  (1983)  closure  model  tc  middle  and  high-level  cloud  systems.  As 
a  first  step  Flatau  (1985)  attempted  to  develop  a  level  2  or  level  2.5 
version  of  the  model  so  that  the  number  of  predictive  equations  in  a 
mixed-phase  cloud  system  could  be  minimized.  A  paper  entitled 


''Comments  on  Realizability  Conditions  for  the  Mellor  and  Yamada  Level 
2.5  and  3.0  Model's"  by  Flatau  and  Cotton,  has  been  accepted  lor 
publication  (see  Appendix  1)  . 

Since  this  stage  of  research,  the  work  has  shifted  more  to 
refinements  in  the  microphysics  and  radiation  schemes  in  RAMS  to  make  it 
more  suitable  for  simulating  cirrus  clouds.  The  results  of  that  work 
are  described  below. 

3.2  Introduction  and  Background 

Upper  level  clouds  and  large  scale  dynamics 

Clouds  composed  of  ice  crystals  (cirrus)  occur  frequently  ir.  the 
upper  troposphere.  They  are  related  to  a  large  scope  of  atmospheric 
phenomena.  These  phenomena  include  systems  of  convective  origin  such  as 
Mesoscale  Convective  Systems  (Chen,  3984;  Maddox,  1981,  1983;  Wetzel  et 
al . ,  1983),  anvils  of  severe  thunderstorms  in  mid-latitudes  (Keymsfielc, 
1986;  Benetts  and  Ouldridge,  1984),  and  tropical  anvils  (Ackerman  et 
al.  ,  1  987  )  . 

Cirrus  clouds  related  to  large  scale  synoptic-type  situations  are 
well-known  to  occur  in  a  variety  of  situations.  Jet-stream  cirrus  Lave 
beer,  investigated  (Conover,  1960),  ana  the  classical  pre-frontal  ice 
cloud  conditions  were  studied  by  mar.y  (Starr  and  Cox,  1980;  Heymsfield, 
1977)  .  A  summary  of  such  extratropical  cyclone  cloud  systems  is  giver, 
ir.  Cotton  and  Anthes  (1987). 

Morphology 

Not  only  are  high-level  clouds  related  to  different  synoptic  scales 
of  motion,  but  they  exhibit  a  wide  range  of  morphological 
characteristics.  Stratiform  ice  clouds  associated  with  Mesoscale 


Convective  Systems  are  quite  common,  but  they  are  often  inhomogeneous  ir. 
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both  the  vertical  and  horizontal  (Yeh  et  al. .  1986),  with  varying 
degrees  of  embedded  convection.  Such  embedded  convective  cells  are 
probably  most  strikingly  exhibited  in  cirrus  uncinus,  clouds  which  were 
subject  of  several  studies  (Heymsf ield ,  1972;  1975a, b).  The  convective 
nature  of  some  high  level  clouds  makes  description  of  their  bulk 
properties  more  complex,  requiring  modelling  of  smaller  scale  phenomena 
within  an  otherwise  uniform,  stratified  system  (e.g.,  Starr  and  Cox, 
1985a, b).  Idealized  cases  of  stratified  and  uniform  upper  level  clouds 
were  modeled  by  Ramaswamy  and  Detwiler,  1986. 

Mlcrophvslcs 

Variations  in  morphological  characteristics  are  closely  related  to 
the  microphysical  constitution  of  such  clouds.  Cotton  and  Anthes  (1987) 
ir.  their  forthcoming  book  on  cloud  physics  summarized  the  microphysical 
properties  as  found  by  different  authors.  This  is  presented  in  Table  I. 

As  summarized  in  Table  I,  cirrus  clouds  contain  ice  crystals  of 
typical  length  0.5  mm  and  in  concentrations  which  vary  widely  from  104 
to  10°  m  .  The  terminal  velocity  of  the  ice  crystals  is  about  50  cm 
s  ,  which  is  of  the  same  order  of  magnitude  as  the  updraft  speeds.  The 
ice  water  content  of  cirrus  clouds  is  typically  0.2  g  m  . 

Just  how  representative  these  properties  are  of  a  typical  cirrus 
cloud  is  uncertain  as  there  is  much  natural  variability  in  ice  cloud 
structure.  For  example,  it  appears  now  that  ice  clouds  exhibit  quite 
complex  structure  as  far  as  ice  crystal  spectrum:  Platt  (1987) 
investigated  cirrus  crystals  in  the  2.5  pm  to  200  pm  size  range  and 
showed  that  there  is  an  abundance  of  the  smaller  particles.  Heymsfield 
and  Platt  (1984)  attribute  up  to  50%  of  cirrus  extinction  to  very  small 
particles  (Fig.  3.1).  The  increasing  abundance  of  remote  sensing  and  in. 
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Average  particle  spectra  in  different  temperature  ranges  fcr 
the  20-300  crystal  maximum  dimension.  (From  Heymsfield  and 
Platt,  1984). 


r- 


0  nm  ma> 
4)  . 


Fig.  3.1b.  As  in  Fig.  3.1a,  but  for  20-180 
(From  Feyrcsfield  and  Platt,  198 


situ  measur events  suggest  that  sub-visible  ice  elements  may  play  a  mere 
important  role  than  previously  thought  (Radke  and  Hall,  1983).  Earr.es 


Table  I.  Microphysical  characteristics  of  middle  and  high  clouds. 


Prooertv  or  variable 

Value 

Reference 

Concentration  of  ice 

105  _  106 

Eraham  and  P. 

crystals  fm-3] 

Spyers-Durar.  (1967) 

2x105  -  5x105 

Keymsfield  (1975a) 

1.0x104  -  2.5  x  104 

Heymsfield  and 
Knollenberg  (1972) 

6.0x105  -  38.0  x  105 

Rvan  et  al.  (1972) 

2.0xl04  -  8.0  x  104 

Houze  et  al.  (1981) 

Length  of  crystals 

up  to  0.17  mm 

Eraham  and  P. 
Spyers-Durar.  (1967) 

0. 6-1.0  mm 

Heymsfield  and 
Knollenberg  (1972) 

length  0.35-0.9  mm 

Heymsfield  (1975a) 

Terminal  velocity  of 

Typically  50  cm  s~l 

Heymsfield  (1975a) 

crystals 

Max  120  cm  s~l 

Ice  water  content 

0.15-0.25 

Heymsfield  and 

(g  nr3) 

Knollenberg  (1972) 

0.15-0.30 

Heymsfield  (1975a) 

0.20-0.50 

Rosir.ski  et  al.  (  1970! 

Frecipitation  rate 

0.5-0. 7  rnn/h 

Heymsfield  and 
Knollenberg  (1972) 

0.1-1 .0  mm/ hr 

Heymsfield  (1975a) 

Updraft  velocity 

1 .0-1 ,5  m  s-l 

(in  cirrus  uncinus) 

Keymsfield  (1975a) 

2.0-10.0  cm  s-1 

(warm  front  overrunning) 

Keymsfielc  (1977) 

20  cm 

(warm  front  occlusion) 

Keymsfield  (1977) 

25-50  cm  s-1 
(closed  low  aloft) 

Keymsfield  (1977) 

Lifetime  of  individual 

15-25  min 

Heymsfield  (1975a) 

cloud 

(1982)  finds  two  distinct  types  cf  thin  cirrus  at 

mid-  and  upper- 

tropospheric  levels.  The  first  type  is  very  common  and  consists  of 
concentrations  of  103  -  i o4  m3  and  no  particles  larger  than  100  urn.  He 
found  this  type  of  size  distribution  70%  of  the  time  while  flying  ir, 
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"clear"  air  at  cirrus  altitudes.  Figure  3.2  shows  large 
concentrations  of  small  ice  crystals  ir.  the  cirrus  clouds.  Liquid  water 
contents  in  thin  cirrus  can  be  expected  to  be  about  an  order  of 
magnitude  less  than  those  in  Table  I. 

The  crystal  habits  are  of  complex  structure  (Kemysfield  and 
Knollenberg,  1972),  with  columnar,  bullet-like  crystal  and  bullet 
rosettes  dominant.  In  the  small  crystal  range,  Platt  (1987)  notices 
that  "there  are  suprisingly  few  differences  between  individual  flights 
and  between  tropical  and  midlatitude  cirrus,  considering  the  wide  range 
in  locations  and  temperatures  in  which  sampling  took  place."  Ke 
observed  Shimizu  (1963)  columns,  needles,  bullets  and  'rimed'  group 
crystals.  Although  there  is  much  similarity  between  crystals  shapes, 
the  existing  diversity  poses  formidable  difficulty  as  far  as  calculating 
their  single  scattering  radiative  properties  (see  Section  2.2  for  more 
details ) . 

Despite  the  difficulty  involved,  some  microphysical 
parameter izations  of  ice  clouds  have  been  attempted.  For  example, 

Cotton  and  his  col laborators  in  the  series  of  papers  (Stephens,  1979; 
Tripoli  and  Cotton,  1982;  Cotton  et  al .  .  1982;  Cotton  et  al. .  1986;) 
developed  a  parameterization  scheme  for  ice  processes  which  is  in 
routine  use  in  the  RAMS  cloud/mesoscale  model.  Other  groups  use  a 
similar  bulk  microphysics  approach  (Hsie  et  al . ,  1980;  Lin  et  al . , 

1980).  This  subject  is  discussed  in  detail  in  Cotton  and  Anthes  (1987). 
Keymsfield  (personal  communication)  developed  a  technique  for  diagnosing 
the  ice  water  content  and  precipitation  rate  of  stratiform  clouds  from 
radiosonde  data.  Such  a  scheme  is  currently  being  tested  for  use  in 


HCAfi's  Community  Climate  Model. 
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Figure  3.3  is  a  schematic  of  the  complex  microphysical  interactions 
that  come  into  play  between  different  water  classes.  This  interaction 
diagram  is  from  Cotton  and  Anthes  (1987)  and  shows  how  aggregates, 
graupel,  ice,  cloud  droplets,  and  rain  droplets  interact  through  various 
physical  processes.  Recently,  Tripoli  et  al.  (1987)  attempted  to 
formulate  a  more  general  microphysical  parameterization  scheme  for  use 
ir.  the  new  RAMS  structure  where  all  these  interactions  could  be 
considered  in  a  consistent  way. 

Radiative  characteristics 

Cirrus  radiative  properties  influence  incoming  solar  and  outgoing 
planetary  radiation  (Freeman  and  Liou,  1979;  Ou  and  Liou,  1984; 

Stephens,  1  980a,b;  Pamar.athan  et  al .  .  1983),  thus  influencing  the  earth 
climate.  It  is  recognized  that  a  change  in  cloud  ccver  will  have  a 
dramatic  effect  or.  the  Earth's  radiative  budget.  Ever,  changes  in 
atmospheric  gases  influence  the  radiative  budget  mostly  through  cloud 
feedback  (Ramanathan  et  al. ,  1983). 

Radiation  plays  an  important  role  ir.  clouds,  and  coupled  with 
turbulence  and  microphysics,  influences  the  dynamics  of  cloud  systems. 
This  was  clearly  established  for  stratocumulus  clouds.  Flatau  (1985) 
and  Chen  and  Cotton  (1983),  as  well  as  other  investigators,  simulated 
off-shore  California  stratus  with  the  coupled  radiative,  turbulence,  and 
microphysical  model,  and  they  have  shown  radiation  influences  on  heat 
and  momentum  fluxes.  In  another  modelling  study  by  Cotton's  group,  Chen 
(1984)  and  Chen  and  Cotton  (1987)  analyzed  the  instability  on  the  top  of 
the  stratiform,  cloud  and  linked  it  with  strong  radiative  forcing.  It 
was  found  that  cloud  top  entrainment  is  triggered  by  the  vertical 
transport  of  turbulence  kinetic  energy  and  is  further  enhanced  by 
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Fig.  3.3.  Schematic  flow  diagram  of  RAMS  microphysics  parameterization 
scheme. 
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evaporative  cooling.  The  energy  driving  the  entrainment  was  initially 
produced  by  cloud-top  destabilization  caused  by  long-wave  radiative 
cooling.  In  addition,  it  was  found  that  the  competition  between  shear 
and  radiative  coding  is  the  key  to  determining  whether  the  cloud  layer 
is  shear-driver,  or  buoyancy-driven.  For  the  case  with  strong  shear,  the 
cloud  top  radiative  cooling  was  found  to  be  balanced  by  entrainment 
heating.  Furthermore,  the  effect  of  shear  may  become  important  when  the 
cloud  top  radiative  cooling  is  weakened  by  large  scale  subsidence,  the 
presence  of  high  level  clouds,  or  a  moist  layer  above  the  capping 
inversion. 

Fadiation  was  also  linked  with  strong  dynamical  effects  in 
Hesoscale  Convective  Systems  through  interaction  with  its  extensive, 
stratiform  anvil  (Chen,  1986;  Chen  and  Cotton,  1987).  Recently,  Tripoli 
(1986;  Tripoli  and  Cotton,  1987)  has  proposed  that  cloud-top  radiative 
cooling  traps  gravity  wave  energy  beneath  the  cirrus  canopy  of  ar.  MCS. 
The  resultant  trapped  waves  play  a  crucial  rcle  ir.  the  up-scale 
development  of  a  KCS  into  a  mesoscale  convective  complex. 

Clouds  also  provide  a  large  atmospheric  target  for  backscattering 
to  lidar  systems  (Radke  and  Hall,  1983)  and  will  therefore  generally 
play  a  beneficial  rcle  ir.  the  operation  of  wind  sensing  lidars  unless 
the  signal  is  too  strongly  attenuated. 

Numerical  modelling 

(a)  Radiation  models 

Up  to  now,  most  (but  not  all)  of  the  cirrus  cloud  models 
concentrated  on  the  radiative  part  of  the  problem  without  much 


consideration  to  dynamics  and  microphysics.  For  example,  Stephens 
( 1980a, b)  used  a  doubling  and  adding  model  to  solve  the  Radiative 


Transfer  Equation  (RTE)  and  presented  results  of  radiative  properties  of 
cirrus  clouds  in  the  IR  region.  He  considered  two  types  cf  ice 
particles,  cylinders  and  spheres,  but  assumed  a  mono-disperse  spectrum 
for  the  broad-band  calculations.  His  results  indicate  that  reflectance 
of  cirrus  clouds  is  small,  but  strongly  influence  cloud  heating 
profiles.  The  reflectance  itself  is  governed  by  cloud  microphysics. 

Roewe  (1977)  and  Roewe  and  Lieu  (1978)  used  the  discr ete-ordir.ate 
method  (DOM/  to  solve  the  RTE  in  cirrus  clouds.  In  their  study,  the 
impact  of  high  level  cirrus  clouds  upon  the  infrared  cooling  of  the 
atmosphere  was  investigated.  Eand  ty  band  cooling  was  calculated  for 
all  major  gaseous  absorption  bands  including  water  continuum.  In 
addition,  the  influence  of  a  low  level  cloud  on  the  infrared  cooling  was 
investigated  in  the  presence  of  cirrus  clouds  of  varying  thickness.  DOM 
is  applied  to  non-iscthermal ,  inhomogeneous  cloudy  atmospheres.  They 
assumed  that  cirrus  clouds  contain  a  mono-disperse  spectrum  of  randomly 
oriented  crystals  in  a  horizontal  plane.  The  crystals  are  long  circular 
cylinders  with  a  mean  length  of  200  ^m,  30  pm  in  radius,  and  0.05  cc'- 
rr.ean  concentration. 

Liou  (1972)  used  a  simple  twe-stream  model  to  calculate  reflecticn, 
transmission,  and  absorption  cf  several  cirrus  cloud  layers  as  a 
function  of  the  solar  zenith  angle  ir  the  visible  and  near  infrarec 
portion  of  the  solar  spectrum.  He  used  single  scattering  properties  cf 
ice  as  determined  by  his  earlier  works  (Liou,,  1  972a, t). 

Welch  et  al .  (1980a, b)  used  a  two-stream  model  to  calculate  the 
variability  of  ice  cloud  radiative  properties  at  selected  solar 
wavelengths.  The  modified  gamma  distribution 
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was  used  for  crystal  size  distribution.  Several  crystal  shapes  were 
assumed,  but  the  treatment  of  non-sphericity  was  rather  approximate  in 
the  calculations.  One  of  the  conclusions  was  that  the  particle  size 
distribution  may  strongly  influence  reflectance,  transmittance,  and 
absorptar.ee  for  wavelengths  between  1.3  pm  and  2.0  pm.  They  suggest 
that  measurements  of  reflectance  ratios  may  provide  information  about 
size  spectra.  Welch  et  al.  (1980b)  notice  that  radiative 
characteristics  are  extremely  sensitive  to  the  asymmetry  factor  g. 

Ir.  addition  to  these  detailed  radiative  studies,  high  level  clcuds 
were  studied  in  relation  to  climate  (Schneider  et  al . ,  1978;  Stephens 
and  Webster,  1979;  Ramanathan  et  al. .  1983). 

(t)  Coupled  Radiative-Dynamical  model 

To  our  knowledge,  there  are  two  recent  models  where  upper  level 
cloud  dynamics  and  radiation  arc  considered.  These  are  described  in 
papers  by  Starr  and  Cox  (1985a, b)  and  Ramaswamy  and  Detwiler  (1986). 
Ramaswamy  and  Detwiler  (1986)  used  a  ID  model  of  monodispersed ,  randomly 
oriented  ice  crystals.  Transfer  of  radiation  is  modeled  with  the  two 
stream  model.  Fifty-six  spectral  intervals  are  considered.  The 
equivalent  sphere  assumption  was-  made  to  obtain  single  scattering 
albedo.  The  cloud  is  located  between  10  and  11  km  in  a  model 
atmosphere.  Simplified  thermodynamics  was  considered,  with  evaporation 
and  condensation  processes  influencing  crystal  evolution.  Cne  of  the 
findings  was  that  radiative  properties  are  dependent  cn  crystal 
concentrations  and  crystal  length.  The  simulated  cloud  resembles 
cirrostratus,  but  the  authors  note  "comprehensive  understanding  of 


cirrus  cloud  will  require  detailed  treatments  of  both  nicrophysics  and 
convection.  " 

A  two-part  paper  by  Starr  and  Cox  (1985a, b)  is  probably  the  most 
comprehensive  to  date.  The  model  has  four  prognostic  equations:  1) 
potential  temperature,  2)  water  vapor  specific  humidity,  3)  ice  water 
content,  and  4)  vorticity.  Phase  changes  are  evaluated  every  five 
minutes.  Bulk  microphysics  parameterization  is  used  with  only  one  class 
of  ice.  Depending  on  ice  water  content  this  ice  class  changes  its 
characteristics  (terminal  fall  speed,  etc.).  Radiation  is  represented 
by  broadband  infrared  and  broadband  shortwave.  The  environmental 
vertical  velocity  can  be  imposed  ar.d  different  environmental  conditions 
(thermodynamics)  are  possible.  The  diurnal  cycle  is  studied.  This  2D 
model  with  a  100  m  grid  predicts  generating  zones  for  cirrus,  and 
dynamics  of  generating  zones  car.  be  quantitatively  studied.  Crystal 
habit  and  size  distribution  are  found  to  be  quite  important  in 
determining  the  overall  water  budget. 

3.3  Summary  of  Work  Accomplished  Usir.g  RAMS 

In  this  section  we  summarize  work  accomplished  using  RAMS  and 
supported  by  the  AFGL.  Two  modules,  radiation  and  microphysics,  and  the 
work  related  to  their  "fine  tuning",  are  described  here.  First  we 
describe  the  current  radiation  module  (Section  3.3.1).  Next  we 
summarize  current  work  on  establishing  more  reliable  reflectance  and 
transmittance  information  for  upper  level  clouds  from  knowledge  of 
single-scattering  albedo  (Section  3.3.2).  Subsequently  new  developments 
into  multiscattering  code  are  described  (Section  3.3.3).  Finally,  the 
generalized  microphysics  package  (Section  3.3.4). 

3.3.1  Current  radiation  scheme 


Details  of  the  present  radiation  employed  ir.  RAMS  has  been 
described  elsewhere  (Chen  and  Cotton,  1983).  We  present  only  a  brief 
outline  of  that  scheme  with  particular  emphasis  on  how  radiation-cloud 
interactions  are  treated  and  on  the  improvements  1  squired  to  apply  the 
parameterization  to  ice  crystals. 

The  current  radiation  scheme  is  based  on  the  adding  method  for 
shortwave  and  the  ercissivity  approach  for  longwave  radiation.  Both 
molecular  gases  and  water  absorption  are  considered.  The  reflectance 
and  transmittance  parameterization  is  based  on  the  two-stream  solution 
v.ith  empirical  modi  f  ica  tior.s  (Stephens,  1  97  8b). 

(a)  Shortwave  parameterization 

The  shortwave  radiation  parameterization  is  based  or.  the  adding 
method.  Only  the  r.et  upward  F4  and  downward  F*  components  of  radiation 
are  considered,  ir. dependent  of  wavelengths.  Fig.  3.4  illustrates  a 
simple  layered  two-stream  model  hared  on  a  model  developed  by  Stephens 
and  Webster  1  1 97 9 )  and  Stephens  (  1978a, b).  For  such  a  model,  the  fluxes 
can  be  written  as 

F*(n  +  1)  =  Rfc(!,n+1)  F4(n+1)  +  V*(n~l/2)  (3.1) 


and 


F’""  ■ 


+  V  t  ( n+l  /  2) 


(3.2) 


where  n  represents  the  layer  under  consideration .  The  convention  that 
is  used  has  r.  increasing  downward  with  greater  optica]  depth.  Fe(r.), 
A(r.)  and  Tr(r. )  represent,  respectively,  the  reflectance,  absorptance  and 
transmittance  in  the  nth  layer.  The  reflectance  ne(l,n-l)  represents 
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the  combir.ed  multiple  reflectance  contributed  by  all  layers  above  the 
(n  +  l)th  layer  and  the  additive  terms  Vn(n+l/2)  represent  the 
contributions  of  sources  of  radiation  multiply  scattered  by  each  layer. 
Details  of  the  algorithms  that  generate  these  multiple  scattering 
factors  can  be  found  in  Stephens  and  Webster  (1979). 

The  major  aspect  of  the  shortwave  parameterization  is  tne 
specification  of  the  individual  layer  reflectances  and  transmittar.ces . 
This  is  a  fairly  straightforward  task  for  the  clear  sky  (e.g, ,  Stephens 
and  Webster,  1979),  but  is  more  complicated  for  a  cloud  layer.  The 
approached  used  for  cloud  is  based  on  that  of  Stephens  (1978b)  who 
employed  a  twc  stream  model  as  a  basis  for  determining  Re(n)  and  Tr(r.). 

The  most  important  parameter  ir.  such  a  model  is  the  optical 
thickness  r 0f  the  cloud  layer.  Estimates  of  range  from  5  <  < 

500.  Formally,  the  optical  thickness  varies  with  wavelength.  Thus, 

<r 

f  f(r)  C(x)nr2drdz 

o  (3.3) 

where  x  -  2rr/X,  f(r)  represents  the  droplet  spectral  density  function 
ar.d  C(x)  is  an  efficiency  factor  for  extinction  which  car.  be  determined 
from  Hie  theory.  Stephens  concluded  from  tis  model  calculations  that 
the  primary  factor  controlling  t,;  is  the  liquid  water  content  (LWC)  of 
the  cloud.  Changes  ir.  the  drop-size  distribution  were  masked  by  the 
large  variations  in  LWC.  He  thus  approximated  ty  as 

tn  ~3 / 2  W/re  ,  (3.4) 

where  rg  (micrometers)  is  an  effective  radius  of  the  droplet 


■c.,(x)  = 


Az 

/ 

0 


distribution  defined  as 


r 

e 


<r 

/  r  Krldr 


aj 

/  r2f(r)dr 


o 


(3.5) 


and  W  is  the  liquid  water  path  LWP  (gm  2).  The  LWP  is  defined  as 


LWP  =  w  =  /  P0rc  dz 


(3.6) 


where  rc  is  the  cloud  water  mixing  ratio.  Two  other  optical  properties 
(in  addition  to  optical  thickness)  are  required  tc  complete  the 

specification  of  Pe  ar.d  Tr.  They  are  the  single  scatter  albedo  cjq  and 

backscatter  factor.  These  are  also  related  to  cloud  liquid  water 
content  in  the  scheme  of  Stephens  (1978b).  Therefore  the  shortwave 
characteristics  of  the  cloud  layer  are  completely  specified  by  cloud  LWC 
and  solar  geometry. 

(b)  Longwave  scheme 

The  longwave  scheme  is  based  on  the  emissivity  approach.  This 
method  is  somewhat  traditional  for  clear  sky  calculations  and  the 
approach  used  follows  that  of  Rodgers  (1967)  with  the  modification 
suggested  by  Stephens  and  Webster  (1  97  9)  to  ac.count  for  the  effects  of 
water-vapor-pressure-dependent  absorption  in  the  spectral  region  8(jcr.  <  X 
<  13  jim.  Clouds  in  the  longwave  region  are  treated  in  the  same  manner 
using  the  concept  of  effective  emissivity.  This  quantity  treats  the 
scattering  process  as  if  it  were  an  absorbing  process  thus  providing  a 


modification  to  the  true  emittance  of  the  cloud.  Using  detailed 
spectral  model  calculations  of  IF  transfer  through  low-level  clouds, 
Stephens  (1978)  parameterized  the  effective  emissivity  as 

eti  =  1  -  exp  (-aQuW)  ,  (3  ?) 

where  W  is  the  liquid  water  path  W  =  /  Pc  r^  dz.  and  r^  represents  the 

mixing  ratio  of  all  condensed  water.  He  evaluated  the  coefficients  for 
upward  and  downward  emissivity  as 

aQ«  =-  0 . 1 3m2g~1  and  aQi  =  0.158m2g_1  • 

The  coefficients  act  and  aQi  correspond  to  a  mass  absorption  coefficient 
which  may  be  compared  to  the  values  of  aT  used  by  Oliver  et  al.  (1978). 
As  mentioned  previously,  Stephens  found  little  variation  ir.  a^t  and  ac^ 
with  changes  in  drop  size  distributions. 

The  combined  long-  and  short-wave  radiative  properties  of  clouds 
are  thus  jointly  derived  in  this  scheme  from  cloud  LWP.  Therefore  both 
the  long-  and  short-wave  radiative  properties  change  in  a  related  manner 
when  LWP  is  varied.  This  coupling  of  the  cloud  radiative  properties 
through  the  LWP  is  an  important  feature  of  the  scheme  enabling  the 
proper  evaluation  of  the  total  radiative  budget  of  any  given  cloud  layer 
and  the  variation  of  this  energy  budget  with  changing  cloud  properties. 

There  are  a  number  of  complications  and  uncertainties  that  arise 
when  applying  this  scheme  to  ice  clouds.  These  are: 

(1)  The  relationship  of  optical  thickness  and  u  to  ice  water  path 

o 

(IWP)  has  not  been  established  and  simple  application  of  these 
relationships  based  on  liquid  water  cloud  is  not  appropriate. 


(2)  The  coefficients  appearing  in  the  emissivity  parameterization  (3.7) 
are  uncertain  for  ice  clouds.  This  uncertainty  is  aggravated  by 
non-negligible  effects  of  scattering  by  ice  particles  in  the  IF 
region. 

(3)  It  is  not  possible  to  use  more  complex  theory  unambiguously  to 
provide  the  required  parameterie  relationships  as  was  the  case  for 
water  clouds.  The  intrinsic  scattering  properties  of  ice  crystals 
in  the  atmosphere  remains  uncertain  due  to  complexities  ir.  shape, 
orientation,  density,  size  and  concentrations.  We  thus  need  to 
rely  more  heavily  on  observations  in  cirrus  and  ice  clouds  (such  as 
will  be  provided  by  the  FIRE  data  set)  to  establish  confidence  in 
the  scheme. 

3.3.2  Single  scattering  for  hexagonal  ice  crystals 
In  the  previous  section  we  described  the  radiative  code  currently 
used  in  RAMS  and  the  aspects  of  the  scheme  that  need  to  be  re-addressed . 
The  vast  majority  of  scattering  particles  in  the  upper  level  clouds  are 
not  spherical,  most  of  the  current  radiation  calculations  for  ice  clouds 
are  based  on  the  assumption  that  the  ice  particle  is  a  homogeneous 
sphere.  This  is  not  realistic  as  a  wide  variety  of  shapes  and  sizes 
have  been  observed  (Heymsfield  and  Knollenberg ,  1972).  Plate  crystals 
and  columnar  crystals  of  hexagonal  cross  sections  are  quite  common.  Ir. 
addition,  the  orientation  is  not  necessarily  random,  with  preferred 
directions  of  falling  governed  by  the  hydrodynamics  of  the  particle. 

Some  account  of  these  complexities  need  consideration  in  developing 
the  radiation  parameterization.  In  this  section  we  discuss  how  this 
problem  can  be  approached  using  anomalous  diffraction  theory.  A  brief 
account  of  this  approach  with  the  new  developments  that  we  have  made 


will  now  be  presented.  The  calculations  discussed  here  have  direct 
application  to  use  in  more  sophisticated  radiation  models  (such  as 
discussed  below)  which  are  to  be  tested  using  observational  data  and 
which  provide  the  ultimate  framework  for  the  parameterization. 

Until  now  a  rigorous  electromagnetic  wave  solution  has  not  been 
obtained  for  hexagonal  cylinders.  Approximate  methods  are  based  or. 
geometrical  optics  (Cai  and  Liou,  1982;  Jacobowitz,  1971),  semi- 
empirical  techniques  (Senior  and  Weil,  1977),  and  approximate  and 
expensive  Extended  Boundary  Condition  Method  (EP^M)  (Mugnai  and 
Wisoombe,  1  980). 

The  alternative  but  still  approximate  theory  is  Dased  or.  the 
principles  of  geometrical  optics  and  diffraction.  This  anomalous 
diffraction  theory  (ACT)  (e.g.,  Stephens,  1984)  holds  for  particles  with 
a  refractive  index  close  to  unity 

(m  -  1)  <<  1 

and  ratio  of  particle  size  and  wavelength  a/X  >>  1  .  Since  the  index  of 
refraction  is  close  tc  1,  the  problem  of  total  scattering  ir.  the  ACT 
theory  reduces  to  calculations  of  the  interference  between  the  almost 
straight  transmission  and  the  light  diffracted  according  to  Huygens' 
principle . 

Thus,  in  the  anomalous  diffraction  approximation,  the  forward 
scattering  amplitude  S(0)  is  given  by  (Bryant  and  Latimer  1969) 

2  • 

S(0)  =  f  (l  -  e_i*  )  dA 

A  (3.8) 


where  9*  =  kd(m  -  1),  d  is  the  particle  thickness,  and  m  is  the  complex 
refraction:  m  =  n  -  in’ .  The  quantity  9*  is  the  complex  phase  shift  of 


light  passing  through  the  particle  relative  to  that  passing  around  it. 


The  kdn'  term  contributes  to  absorptive  attenuation. 

The  extinction  and  absorption  coefficients  are  defined  as 


Cext  =  I  Re  U  (  1  ‘  ** 


)  dA 


(3.9) 


and 


r  _  I  r  , ,  -2kdn' .  . . 

abs  =  A  (  {1  "  e  >  dA 

A 


(3.10) 


The  single  scattering  albedo 

_  Cext  _ _ ^abs 
w  C 

abs  (3.11) 

represents  the  probability  that,  if  a  photon  interacts  with  an  elemental 
volume,  it  will  be  scattered  rather  than  truly  absorbed. 

Here  we  present  calculations  of  ^ex^.,  C^g,  and  w  for  hexagonal 

crystals . 

The  geometry  of  the  problem  is  illustrated  in  Fig.  3.5. 

The  sun’s  rays  are  coming  perpendicular  to  the  2-3  face  (=  -  90°) 
in  case  1  (Fig.  3.6)  and  perpendicular  to  face  1-2  ( «  =  60°)  in  case  2. 
Only  two  cases  are  considered,  but  they  represent  extreme  situations. 

The  length  of  crystal  is  c  and  the  length  of  distance  1-2  is  a. 

Scattering  geometry  for  Case  1  is  presented  in  Fig.  3.6.  We  have 


d(y)  =  2\ |  3y  on  dQdj 


(3.12) 


Extinction  efficiency  is  given  by 


Cext  =  t  2c  /  f<y)  dy  ’ 
d0 


(3.13) 


where  f(y)  =  1-e  ° ^ cosk(n-l )d (y)  .  The  factor  2  appears  because 

of  symmetry,  and  the  real  part  of  (2)  is  taken,  A  =  2ac.  Length  a(y)  is 
one  when  the  sun  ray  travels  inside  the  crystal. 

For  case  2  geometry  we  have  on  d^^ 

d(y)  =  a  +  —2—  y 

\  |  3  (3.14) 


and  extinction  efficiency  is  given  by 


,  1 

Cext  =  A  2C  /  f(y)  dy  • 

°0  (3.15) 

After  some  algebra,  it  can  be  shown  that  the  extinction  efficiency 
factor  for  case  1  is 


ext 


=  2 


_ 2tan& _  e.  ptanP  _  g-p tang cos(p) tan2 P 

p\ |  3sec^p  sec2g  sec^g 


-l/2ptang\|  3 


'.ft 


The  absorption  coefficient  (case  1)  is 


Q  _  l^ptanp  _  1 _ 1 _ +  ig-ptanpv  |  3 _ 1 _ 

abs  ~  2  2ptar,p  \  |~3  2  ptanp\|T  (3.17) 


Similarly,  for  case  2,  we  have 


.  -l/2ptanp  -l/2Ptanp 

Cext  =  2\  - p - Sin2p  -  - cosptanf 


.  -ptanp  .  -ptanp  , 

- sinp  +  “ - cosptanp  +  tan^p  +  1 


/  (tan  P+1)  , 

(3.38) 


and 

p-ptar.p  fi-2ptar.p 

^“abs  “  ~  ptar.p  ptanp  2  (3.19) 

The  results  are  shown  in  Figs.  3.7-3.10.  In  Fig.  3.7,  0  (curves 
which  for  large  values  of  p  converge  or  oscillate  around  2)  and  are 

shown  for  different  values  of  complex  parts  of  the  refractive  index.  Ir. 
Fig.  3.8  we  compare  extinction  for  both  cases  and  Fig.  3.9  shows 
averaged  coefficients.  Fig.  3.10  shows  single  scattering  albedo  as  a 
function  of  phase  shift  for  different  values  of  refractive  index. 

3.3.3  Multiple  scattering  code  -  matrix-exponential  method 
We  now  outline  the  sophisticated  method  of  solving  the  radiative 
transfer  equation.  This  method  has  many  advantages  over  other 
techniques.  On  the  one  hand,  it  will  be  employed  in  its  full 
sophistication  to  test  the  appropriateness  of  the  scattering  theories 
against  observations.  On  the  other  hand  and  because  of  the  analytic 
nature  of  the  method,  the  theory  can  be  simplified  down  to  the  level  of 
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of  Case  1  and  Case  2  extinction  coefficient  for 
5. 


SI 

a  two  streair.  model  as  used  in  the  original  RAMS  scheme.  We  now  outline 
the  method. 

Introduction 

Recently  Flatau  and  Stephens  (1987)  reviewed  and  extended  the 
matrix-exponential  technique  for  solution  of  the  radiative  transfer 
equation  (Aronson,  1  972;  Karp  et  al .  .  1980;  Waterman,  1981).  The 
matrix-exponential  technique  allows  one  to  write  the  solution  of 
radiative  transfer  equations  for  the  homogeneous  layer  as 

T(x)  =  6At  I (o )  +  Inh  (3.20) 


(3.22) 


but  the  technique  is  very  general.  Two-point  boundary  conditions  are 
considered  using,  as  an  example,  the  two-stream,  one-level  model.  Ar, 
algorithm  is  proposed  to  solve  RTE  in  terms  of  matrix  multiplication 


A1  Atj  A2  At2  AnArn 

i  e  . . .  .e 


(3.23) 


Such  a  method  avoids  finding  "continuity”  constants,  thus  saving  time 
and  computer  memory. 

We  show  that  the  comutator 


[  A ,  R]  =  AB  -  BA 


(3  .24) 


plays  ar.  important  role  in  RT  and  that  it  expresses  the  fact  that  the 
two  layers  interact  (add)  in  a  non-linear  way.  The  doubling  technique 
(Grant,  1969a, b,c)  is  also  discussed,  and  the  identity 


A2tN  ,  AtN. 2 
e  =  (e  ) 


(3.25) 


is  used  to  formulate  the  doubling  scheme. 

The  star  product  is  shown  to  be  equivalent  to  standard  matrix 
multiplication  of  the  form 


AjTj  A2t2 
e  e 


(3.26) 


We  also  hypothesize  about  the  role  of  matrix  operators 


exp (/  A(r)dt) 


and 


S3 


exp (A  /  w- (r)dt) 

o  (3.27) 


for  inhomogeneous  layers.  It  is  suggested  that  they  may  play  an 
important  role  ir.  the  approximate  solution  of  RTE  or  PTE’s  theoretical 
studies . 

Instead  we  formulate  the  problem  in  the  short  sub-section  below. 

Matrix-exponential  scheme  -  formulation  of  the  problem 
The  equation  describing  the  transfer  of  monochromatic  radiation 
through  a  plane  parallel  homogeneous  atmosphere  is  given  by  (Stamnes, 

1  986  ) 

M  -  J(t,m,0).  (3  2 


where  Kt.p.p)  is  the  specific  intensity  at  level  t  ir.  a  cone  of  unit 
solid  angle  along  the  direction  p.9,  T  is  the  optical  depth,  <f  the 
azimuthal  angle,  p  the  cosine  of  the  polar  angle,  and  J  is  a  source 
term. 

The  discrete  ordinate  approximation  to  the  zero-order  Fourier 
component  of  (3,28)  can  be  written  as  (Stamr.es  and  Swanson,  1981) 


M 


i 


dli(t 

dt 


Ii  (t) 


j,  4  cijV  -  °i(0' 


(3.29) 


where  I^t)  =  IU.pJ,  Qi(r)  =  C(T,pi>,  and 


(3.3  0) 


2n-l 

cij  =  l/2w0  aj(21  +  l)g1P1(Mi)P1(Mj)  . 


The  system  (3.29)  of  2n  coupled  differential  equations  car.  be 
written  in  matrix  form  as 


+1 

dl_ 

-f.  * 

dr 

-a  -P' 

I 

dll 

Pa 

l“ 

dx 


I (t ) 
I(x,±p2) 


I(x,+Jln) 


where  the  elements  of  the  r.  x  r  matrices  a  and  p  are  given  by 


a 


ij 


-  6  .  ;  ) 
1  J 


6  .  . ) 
1J 


(3.31) 


=  1  if  i  =  j  and  =  0  if  i  4  j. 

-1- J 

We  wish  to  solve  the  radiative  transfer  problem  for  a  plane  parallel, 
planetary  atmosphere  having  all  its  inhomogeneities  confined  to  the 
vertical  direction.  Figure  3 .4  illustrates  the  problem  and  introduces 
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some  of  the  notation  conventions  to  be  used.  Each  layer  in  the 
atmosphere  is  assumed  to  be  homogeneous  but  may  have  an  arbitrary 
optical  thickness.  At^  ~  ~  Ti-l*  Eoal  to  compute  the  angle- 

dependent,  specific  intensity  at  all  points  in  the  atmosphere. 

I  <r;  4,9),  where  n  is  the  consine  of  the  zenith  angle,  and  9  is  the 
azimuth  angle  measured  from  the  sun's  meridian. 

Ke  may  write  the  discrete  ordinate  equations  (2.23)  as 


dl  = 
dt 


AI  +  S 


(3.32) 


where  S  is  the  source  vector.  The  form  of  the  matrices  appearing  in 
these  equations  is  not  important  for  the  current  discussion.  Writing 
Eq.  (3.3  2)  ir.  the  form 


dl 

dt 


AI  =  S 


and  multiplying  by  e 


-At 


we  obtain 


ArI) 


-  At^ 


(3.33) 


Integrating  between  0  and  t,  the  above  becomes 


e  AtT(t)  =  1(0)  +  f  e_ArS( r)dt 


(3.34) 


Finally , 
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T(t) 


eArI( 0)  +  /eA(T-t) 
o 


S(t )dt 


(3.35) 


Here  we  assumed  that  A  is  constant.  In  more  specific  form  the  general 
solution  (3.35)  can  be  written  for  one  homogeneous  layer  as 


I(xA)  =  e  Ht.^)  +  if- 


where 


ynh 

i 


Ti 

/ 


Ai(ti-1 


T 


i-1 


t)  _ 

£(t )  dt 


(3.36) 


If  we  are  not  interested  in  the  interior  radiation,  the  system  of  Eq. 
(3.35)  can  be  simplified  immediately  to 


:<t=xn)  = 


ANAtN 


AN-1Ai:N-1 


.e 


AjAt  j 


1(0) 


+  7* 


From  the  solution  (3.35)  we  can  form  the  matrix  equation 


I(t)  «  QI(0)  +  I1^ 


where 


Q  =  eAt. 


(3.37) 


Using  definition  (3.37)  of  transfer  array  Q  we  can  introduce  G-matrix 


for  the  entire  atmosphere  as 
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°2  Q1 


(3.38) 
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► 
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i 
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where 


Q. 

l 


AiAti 

e 


(3.39) 


A  most  aggravating  property  of  the  matrix  exponential  is  that  the 
familiar  additive  law  fails  unless  we  have  commutativity  (Moler  and  van 
Loan,  1978) ,  i.e. , 

etAeTp  =  et(A+B) 

if  and  only  if 


It  is  evident  that  if 


[A,R]  -  AB  -  BA  =  0. 


t  C. 


Q.)  = 


(3.40) 


for  every  combination  of  i,j  then  Q  in  Eq.  (3.38)  could  be  expressed  by 


/  NA  dt 
C  =  e  c 


(3.41) 


Eecause  (3.40)  is  not  satisfied  it  is  equally  evident  that  commutators 
will  play  an  important  role  in  the  solution  to  RTE.  '  Indeed  it  is  known 
that  adding  slab  reflectances  and  transmittances  does  not  follow  linear 
law,  which  gives  a  physical  meaning  to  non-commutativity. 

It  should  be  emphasized  that  equation  (3.37)  is  only  a  formal 
solution,  since  we  are  dealing  with  a  two-point  boundary  value  problem 


iR 


in  which  T~ (0)  and  I+(t„)  are  to  be  determined.  However,  we  can 
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overcome  this  problem  and  show  equivalence  of  the  matrix-exponential 
technique  with  the  doubling  and  adding  scheme. 

3.3.4  Unified  approach  to  microphysics  scheme 

Radiation  properties  as  described  in  Section  3.3.2  are  related  to 
single  crystals  only.  To  use  it  in  the  RAMS  model  we  have  to  weight  the 
single  scattering  albedo  with  the  crystal  size  distribution  to  obtain 
the  bulk  radiation  characteristics.  Currently  the  microphysics  scheme 
is  based  on  the  Marshall-Palmer  distribution,  but  we  expect  that  some 
other  distributions  may  fit  the  experimental  results  better.  For 
example  recent  experimental  results  (FIRE  data)  as  well  as  new  analysis 
by  Heymsf ield  and  Platt  (1984)  indicate  that  small  crystals  may 
contribute  up  to  50%  of  overall  radiative  characteristics  of  ice  clouds. 
It  strorgly  indicates  the  necessity  of  introducing  a  twc-motial 
distribution  function  with  the  primary  mode  describing  small  ice 
crystals,  numerous  in  concentration  but  small  in  water  content,  and  a 
secondary  mode,  containing  large  crystals  with  small  concentrations  tut 
large  ice  water  content.  Below  we  summarize  a  new,  unified  approach  to 
the  bulk  microphysics  scheme.  We  give  examples  how  to  efficiently 
calculate  the  generalized  moments  for  four  different  distributions.  It 
is  then  easy  to  obtain  a  formal  solution  for  the  bi-modal 
parameterization. 

This  summary  is  an  outgrowth  of  several  previous  works  by  Cotton 
and  his  co-workers.  These  include  rain  parameterization  from  Mantor.  and 
Cotton  (1977)  and  Tripoli  and  Cotton  (1980);  ice  parameterization  from 
Cotton  et  al .  (1982)  and  Stephens  (1979);  and  aggregation 


parameterization  from  Cotton  et  al.  (1986). 


Several  suggestions  as  to  the  type  of  distribution  to  be  expected 
in  the  clouds  are  considered  here.  The  Marshall- Palmer  function  is  used 
because  it  is  a  basis  of  current  implementation.  AFGL's  truncated 
distribution  (Eerthel  and  Plank,  1983)  is  discussed  as  well  as  Clark  and 
Hall’s  (1983)  log-normal  function.  The  gamma-type  function  is  also 
considered . 

Size  spectra 

Several  empirical  formulas  have  been  proposed  for  the  spectra  of 
various  water  categories.  The  spectral  density  function  r.(D)  is  defined 
such  that  r.(E)dD  is  the  total  number  concentration  of  particles,  whose 
diameter  is  between  C  and  D+dD,  per  unit  volume  of  air.  Hence, 

D 

N(P)  =  /  n(D)dD, 

0  (3.42) 

is  the  total  concentration  of  particles  of  diameters  smaller  than  D,  per 
unit  volume  of  air. 

(a)  Gamma  distribution 

A  family  of  gamma  distributions  was  considered  by  several  authors 
(see  e.g.  Scott  (1968)).  The  defining  equation  for  this  family  is 


\H1 


N(D) 


r  <\m)d 


<f>y 

m  m 


exp (-  j^) 
m 


(3.43) 


where  D  is  the  diameter  of  particle,  N  the  j_g  total  concentration,  P_ 

t  m 

is  mean  diameter,  a  and  \)  are  coefficients  and  is  the  gamma  function 
(b)  Marshall- Palmer  distribution 
The  special  case  of  (3.43)  is  often  used  in  cloud  physics. 
Assuming  a  =  1 ,  \)  =  0  we  have 


N  ( D )  =  exp  ( -  jp- ) 


(3  .44) 


The  mean  diameter  is  related  to  the  slope  of  the  distribution  -  (X.) 


(3.45) 


(c)  Doubly-truncated  Marshall-Palmer  distribution 
AFGL  measurements  (Berthel  and  Plank,  1  983;  Berthel  et  al .  .  19  87) 
indicate  that  it  is  often  advantageous  to  approximate  the  spectral  size 
distribution  with  the  Marshall-Palmer  probability  but  truncated  on  both 
sides  of  the  spectrum,  i.e., 

MIDI  -  MP<- 

O  m  (3.46) 


f0r  Dmin  *  ^  *  P'm ^ v *  sr.d  0  otherwise,  where 
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-D  •  / D  -D  /D 

_  min  m  max  m 

S  =  e  -  e 

(d)  Log-normal  distribution 

Clark  (1976)  and  Nickerson  et  a) ■  (1986)  have  advocated  use  of  an 
asymmetrical  distribution  of  the  log-normal  type.  It  is  defined  by 


N(D)  =  - . — exp  [-  (In  D  -  ^)2/(2o2)] 

\|  2jtctD 


(3.47) 


where  \i  and  o2  are  the  mean  and  variance  of  InD,  respectively.  For 


example , 


M 


M  =  In  D  =  lnD 
m 


(3.48) 


where  is  the  geometric  mean  of  the  D's.  We  also  have 

a2  =  (lna  )2  =  ( In if -  In  D  )2 
m  m 


(3.49) 


Moments  and  distribution-weighted  functions 
In  the  development  of  the  bulk  parameterization  scheme,  such  as  we 
are  pursuing  here,  one  is  interested  in  the  average  properties  of  a 
hydrometeor  class.  For  example,  such  a  property  can  be  the  average  fall 
velocity  that  is  "typical"  for  all  size  droplets,  typical  collection 
efficiency,  etc.  Therefore  we  would  like  to  be  able  to  calculate 

averages  f 


f  = 


'max 
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mm 


f  ( D )  N  ( D  )  d  D . 


(3.50) 


Such  integration,  for  an  arbitrary  f(P)  and  K(D),  is  difficult.  For  f 
of  the  form 

f(D  =  KDa ,  (3.51) 

where  K  and  o  are  constants,  and  four  distributions  considered 
previously,  it  can  be  done  in  closed  form. 

If  the  most  general  form  of  f  is  (3.51)  or  the  linear  combination 
of  such  functions,  the  problem  reduces  to  finding  the  generalized  moment 


of  distribution  N(D) 


D 

max 

M°  =  /  ra  N(D)dD. 

D  . 
mir. 


(3.52) 


(a)  Moments  and  averages  of  R-distritution 
From  (3.52)  and  (3.43)  we  have 


(3.53) 


Let  y  =  — ,  then  (3.53)  can  be  written  as 
m 

(3.54) 


From  table  of  integrals 

N  D  » 

m  aaR\Kl) 

(3.55) 


For  f  of  the  specific  fcri 


we  get 


f  =  KDa 


f  =  N. f(D  )F (a)  , 
t  m 


(3.56) 


where 


F(o)  = 

aaR\Kl) 


(3.57) 


(b)  Moments  and  averages  of  Mar shall- Palmer  distribution 
Setting  \)  =  0  and  a  =  1  in  the  expression  for  the  gamma 
distribution  (3.56)  we  get 


M°  =  N.p  “  rta+1) 

t  m 


(3.5  8) 


f  =  Ntf(Dm)F(a). 


(3.59) 


where 


F(a)  =  Ha  +  1) 

(c)  Moments  and  averages  of  doubly-truncated  Marshall 
Palmer  distribution 
From  (3.46)  and  (3.52)  we  have 


max  K. 

M«  =  /  Ca  r-f-  exp  (-  jf 


n  1  SD 

D  .  m 

min 


(3.60) 


Let  y  =  p— .  We  can  write  (3.60)  as 


f,’tDffla  F  /  ya  exp(-  y)dy 
D_ .  _ 


(3.61) 


and  integrating  we  get 


n  m 


The  distribution-averaged  f  is 


(3.62) 


f  =  Ntf(Dm)F(a) 


(3.63) 


where 


F  ( a  ,  D  r  •  ,  T  ) 
max’  mir.  m 


P  D  ■ 

max .  ,  n  ir. , 

y  (a+1  ,  D  )  -  y (a+1  ) 

_  n  ‘'n 


lf  ^Kax'^E  ls  +  ar8e  and  DEin/Dn  is  small ,  further  sirr.pl  if  icaticr.  is 
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F(a  )  =  R  a  +  I ) 


(3.64) 


It  seems  ( Bar.te ,  private  communication)  that  conditions  wh 
snail  are  often  satisfied. 
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(d)  Moments  and  averages  of  leg-normal  distrifcutic 


From  the  definition  (3.52)  we  have 
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Ma  =  /  Da 
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exp  [-( lr.D  -  (i  )2/  ( 2o2 )  ]dD 


(3.65) 


and  let 


v  =  lr.  D.  r-k 

\  \  ~2a 


After  substitution  we  get 
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(3 .66) 


Using  the  expression  for  the  stove  integral  we  have 


M  =■  V 


2,2 

M°  exp(°Y^ 


(3  .67) 


Using  the  above  expression  we  obtain  the  distribution- averaged  function 
of  the  general  form  f  =  K  Da. 


f  =  KtKeM°  exp( 


2  2 
a__a_\ 


(3.68) 


Eut  from  (3.48)  we  have  D  =  e^,  i.e. 

m  c  ’ 


Ke,ia  =  KD  a. 
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Finally,  we  can  express  f  by 


(3.69) 


f  =  Ntf(Dm)F(a.o) , 


where  h'^  is  total  concentration,  f(Dm)  is  a  function  evaluated  for  the 
geometric  mean  of  particle  diameter,  F(a,o)  is  a  function  depending  on 
a  and  variance  a  of  distribution 


F(a,a )  =  exp  (g~^~) 


(3.70) 


Summary  of  distribution  functions  and  their  moments 
We  considered  up  to  now  four  different  distributions:  gamma, 
Marshall-Palmer,  truncated  Marshall-Palmer  and  log-normal.  We  were  able 
to  derive  distribution-averaged  formula  for  the  f  =  KDa,  where  K  ar.d  a 


are  constants  to  the  general  expression 


?  •  (3.7 

where  N  is  total  concentration,  f(D_)  is  the  value  of  the  functicr,  at 
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the  mean  diameter,  and  F(a)  function  of  a  and  possibly  some  other 
parameters.  The  form  of  the  F-function  is  summarized  below. 

GAMMA  DISTRIBUTION 

F(a)  = 

aart\M) 

MARSHALL-FALMER 


F(a)  =  fTa  +  1) 


DOUBLY-TRUNCATED  MARSHALL-PALMER 


It  follows  from  (3.71)  that,  if  f  is  proportional  to  D°  or  if  it  is 
a  linear  combination  of  such  functions,  we  can  formally  perform  the 
calculation  for  one  of  the  distributions  and  substitute  a  different  F(a) 
when  a  different  size  spectra  are  needed. 

This  generalized  approach  was  already  used  to  develop  a  unified 
approach  to  microphysics  parameterization  (Tripoli  et  al . ,  1987). 

In  summary,  the  foundations  of  the  turbulence,  radiation,  and 
microphysical  modules  have  been  developed  to  extend  RAMS  to  the 
simulation  of  cirrus  clouds.  Much,  work  remains  to  complete  this  model, 
however . 

4.0  FORECASTING  AND  DIAGNOSIS  OF  EXTREME  WIND  SHEARS  AND  TURBULENCE 
A  summary  of  work  completed  relating  to  extreme  wind  shears  by 
severe  convective  storms  is  presented  here.  Over  the  course  of  the 
contract  the  focus  has  gradually  shifted  from  the  convective  scale 
downdraft  structure  and  its  influence  on  severe  wind  shear  to  the  upper 
end  of  the  meso-p  (20-200  km)  scale.  It  is  now  apparent  that  long  lived 
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severe  surface  outflow  events  lasting  on  the  order  of  six  to  nine  hours 
may  accompany  a  large  percentage  of  mesoscale  convective  complexes 
(approximately  30  percent  of  the  MCCs  studied  thus  far).  Because  these 
systems  generally  peak  in  intensity  after  midnight  this  represents  a 
genuine  concern  for  the  general  public,  local  forecasters  and  aircraft 
concerns,  such  as  aircraft  tied  down  at  airport  facilities,  in  the 
affected  areas.  A  natural  question  that  arises  then  is  why  do  some 
systems  continue  to  produce  severe  surface  winds  through  the  early 
morning  hours  while  others  peak  in  the  early  evening  and  then  rapidly 
subside?  Also  are  there  signatures  in  the  synoptic  environment  which 
may  allow  one  to  accurately  predict  these  events  ir.  time  to  give  ample 
advice  to  those  who  may  need  the  data  to  either  protect  vulnerable- 
machinery  or  alter  planned  activities? 

Insights  to  the  above  questions  have  been  provided  through  both  the 
numerical  and  observational  studies  conducted  under  the  current 
contract.  A  series  of  numerical  studies  performed  by  Knupp  (1985) 
provided  information  on  the  structure  of  downdrafts  during  the 
developing  and  mature  stages.  He  identified  the  major  downdraft 
circulation  branches  in  Fig.  4.1  as  consisting  of  both  a  mid-level  and 
an  up--down  downdraft  branch.  The  mid-level  branches  are  thought  to 
predominate  in  the  early  stages  of  the  downdraft  when  large  shears  and 
strong  surface  outflows  my  exist.  The  up/down  downdraft  branch  may 
become  a  more  prominent  feature  during  the  mature  storm  stage.  A  factor 
controlling  the  strength  of  this  downdraft  branch  is  the  stability  of 
the  ABL  lying  ahead  of  the  storm.  The  stability  ir  turn  my  be  affected 
by  outflow  layers  from  previous  storm  systems  or  cloud  shadowing  and/or 
evaporation  of  precipitation  falling  into  inflow  layer  ahead  of  the 


Fig.  4.1.  (a)  Schematic  illustrating  primary  relative  flow  branches 

comprising  the  low-level  precipitation-associated  downdraft 
located  along  the  upshear  flank  with  respect  to  the  updraft. 


PRESSURE  (kPb) 


f'j.g.  4.1.  Continued  (b)  Projection  of  primary  downdraft  flow  branch 
onto  a  vertical  east-west  plane.  Physical  processes  along 
each  branch  are  portrayed.  (c)  Sounding  illustrating  the 
relationship  of  downdraft  properties  to  environmental 
structure. 
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storm.  He  stowed  that  downdraft  speeds  along  this  branch  may  attain 
downburst  intensity. 

Other  factors  controlling  downdraft  intensity  were  also  documented. 
The  cloud  model  results  indicate  that  precipitation  size  and  phase 
exerted  a  strong  control  on  the  modeled  downdraft  strengths.  Knupp 
reported  that  melting  accounted  for  approximately  60  percent  of  the 
cooling  in  environments  supporting  low  cloud  bases.  As  cloud  base 
height  increased,  the  relative  cooling  by  melting  decreased  in 
importance.  Sensitivity  to  the  characteristic  size  of  the  precipitation 
was  shown  when  downburst  type  outflows  were  produced  by  merely 
quartering  the  precipitation  size  used  in  the  simulation.  Other 
controls  on  outflow  strength  may  arise  due  to  the  modification  in  the 
AEL  from  either  short- or- long  lived  convective  elements.  Knupp  and 
Jorgensen  (1985)  suggest  that  convective  scale  perturbations  on  a 
resoscale  surface  flew  field  produced  by  a  cor.vectively  induced  mesc- 
high  may  be  cne  mechanism  which  could  sustain  dewnburst  type  outflows 
over  a  long  period  cf  time. 

It  has  long  beer,  recognized  that  certain  radar  echo  configurations 
such  as  bow  echoes  tend  to  be  associated  with  downburst  activity 
(Fujlta,  1  976).  A  Doppler  radar  observation  of  a  bow  echo  system  or.  the 
2nd  of  August  that  occurred  in  Montana  showed  that  the  system  had  both 
the  mic-level  and  up/down  sources  to  the  downdrafts.  This  system  had 
peak  reflectivities  of  70  dBZ  and  peak  outflows  exceeding  30  m/s.  On 
the  storm’s  upshear  flank  at  mid-levels  was  a  well  defined  inflow  toward 
the  storm,  shown  in  Fig.  4.2,  that  was  initially  thought  to  contribute 
significantly  to  the  surface  outflow  through  the  vertical  transport  of 
the  higher  momentum  air  to  the  surface  within  the  predominate  downdraft 


zones.  Later  trajectory  work  indicated  however,  that  the  primary  source 
region  for  the  downdraft  air  was  from  the  downshear  region  of  the  storm. 
The  analysis  indicated  that  the  "rear  inflow  jet"  was  a  mixture  of 
parcels  from  several  different  locations  within  the  storm  including  the 
leading  edge  of  the  system.  Tajectories  released  in  the  rear  inflow  did 
not  reach  the  surface  nor  did  they  approach  the  leading  edge  of  the 
storm.  We  thus  concluded  that  strong  surface  winds  associated  with 
organized  squall  lines  may  not  initially  require  mid-level  transport 
from  the  upshear  portion  of  the  storm  system.  The  analysis  indicatec 
that  the  areal  extent  and  the  strength  of  the  rear  inflow  increased  with 
time  and  may  have  had  a  more  direct  impact  on  the  surface  flow  momentum 
at  a  later  time  in  the  storm's  life  cycle. 

The  surface  analysis  indicated  a  large  mesc-high  of  nearly  6  mb 
behind  the  squall  line.  The  strongest  winds  at  the  surface  were  located 
on  the  southeast  flank  of  the  mesc-high  and  to  the  north  of  the 
strongest  cell  with  the  squall  line-.  The  Doppler  analysis  placed  a 
large  15  r./s  ccwr draft  on  the  northern  flar.k  of  this  cell.  Thus.  ir. 
this  case,  the  strongest  surface  winds  occurred  in  a  region  where 
convective  scale  perturbations  were  superimposed  or.  a  background  flew 
provided  by  the  mesc-high,  supporting  the  hypothesis  given  by  Knupp  and 
Jorgenson  (1985). 

The  storm  of  2  August  was  unusual  in  the  sense  that  it  propagated 
into  the  cold  stable  outflow  wake  of  a  previous  mesc-f)  system  (see  Fig. 
4.3).  Aircraft  measurements  taken  ahead  of  the  squall  line  revealed  a 
elevated  inflow  source  existed  above  the  low-level  stable  layer.  The 
elevated  inflow  layer  in  turn  was  a  remnant  of  the  pre-storm  atmospheric 


boundary  layer  and  was  being  continually  replenished  by  the 


Fig.  4.3.  Surface  mesoscale  features  at  2030  MDT  2  August  1981.  Thir. 

arrows  represent  surface  flow  while  bread  dashed  arrows 
represent  synoptic  flow  above  the  surface  cold  dome  produced 
by  the  storm  labeled  A.  Stippling  denotes  30  aBZ  storm 
contours.  Storm  E  is  the  squall  described  ir.  the  text  and 
the  numbers  are  representative  surface  values. 


synoptic/regional-scale  flow.  This  is  one  mechanism  by  which  storms  may 
maintain  their  intensity  and  could  remain  a  threat  long  after  the  peak 
daily  heating  has  subsided. 

The  strongest  cell  was  located  near  the  vertex  in  the  bow  echo  and 
nearly  80  km  from  the  surface-based  high  B  air  seen  in  Fig.  4.3.  A  more 
complete  analysis  of  the  convective  structure  indicated  that  the 
strongest  cell  was  an  example  of  a  weakly  evolving  supercell  (Foote  and 
Frank,  1983).  The  overall  structure  of  the  squall  line  suggests  that 
the  rotational  characteristics  extended  beyond  the  single  supercell 
within  the  line  to  include  the  ordinary  cells  located  along  the  southern 
50  Ion  of  the  line.  It  is  suggested  that  subclasses  of  convective 
organization  cr  length  scales  of  50  to  100  km  may  exist  along  a  squall 
line  depending  on  the  environmental  conditions  such  as  the  strength  of 
the  wind  shear  and  perhaps  the  resultant  shape  of  the  outflow  produced 
by  different  cells  along  the  line.  These  subclasses  may  be  a 
combination  or  ordinary  and  supercell  type  storms  acting  as  a  coherent 
entity.  It  appears  that  the  Lew  echo  may  represent  one  such  subclass  or 
arrangement  of  convection  along  a  squall  line. 

A  recent  paper  by  Rotur.no  et  al  .  (1987)  indicates  that  a 
controlling  factor  in  the  behavior  of  the  squall  line  is  the  degree  of 
cancellation  between  the  horizontal  vorticity  produced  by  the  solenoid 
associated  with  the  cold  outflow  and  the  horizontal  vorticity  present 
due  to  wind  shear  in  the  inflow  layer.  The  optimal  situation  for  long- 
lived  lines  occurs  when  the  twe  sources  completely ' cancel  each  other. 
Extending  their  analysis  to  a  curved  outflow  pool  one  might  then  expect 
the  optimal  situation  to  exist  only  in  a  relatively  narrow  zone  since 
the  horizontal  vorticity  vector  associated  with  a  circular  cold  pool 


rotates  through  360  degrees  while  the  inflow  horizontal  shear  vector 
retains  its  initial  orientation.  Thus  one  might  expect  a  spectrum  of 
convective  behavior  along  a  curved  outflow  pool  with  the  longest  lived 
features  occurring  in  the  region  where  the  vorticity  vectors  from  each 
source  lie  opposite  to  each  other.  We  speculate  that  this  may  have  hac 
a  control  in  locating  the  supercell  in  the  2  August  case  near  the  vertex 
in  the  bow  echo  where  the  tightest  gradient  was  observed  in  the  circular 
surface  cold  pool  that  accompanied  this  system. 

The  2nd  of  August  squall  line  was  an  integral  component  of  a  KCC 
that  moved  southeastward  over  the  Dakotas  during  the  early  morning  hors 
of  the  3rd.  The  resultant  KCC  was  accompanied  by  severe  wind  reports 
throughout  the  period.  We  have  begun  tc  compile  statistics  on  other 
KCCs  that  also  prcduce  severe  winds  into  the  early  morning  hours. 
Preliminary  findings  are  that  roughly  30  percent  of  the  KCC  episodes 
produce  severe  winds  for  several  hours  after  midnight  with  the  strorgest 
wind  events  exceeding  45  m/s  in  nor- tcrr.adic  activity. 

It  is  our  hcp-e  tc  construct  a  useful  conceptual  model  of  these 
events  based  on  the  work  completec  ir.  part  by  the  funding  provided  tc 
date  from  the  current  contract.  Any  successful  model  will  need  to 
incorporate  the  multi-scale  aspects  of  the  problem  as  touched  or.  above. 
These  include  the  controls;  on  downdraft  strength  uue  to  er.v irormental 
and  microphysical  framework  withir.  which  the  storm  evolves  as  presentee 
by  Knupp  (1985).  Mesoscale  organization  of  large  outflow  pools  and  the 
storms  ability  to  decouple  from  the  surface  may  also  be  important 
factors  in  long  lived  events.  A  better  understanding  of  the  rear  inflow 
jet  and  it’s  ultimate  control  on  the  surface  outflow  is  also  needed. 
Doppler  studies  performed  to  date  indicate  that  a  common  characteristic 


of  bow  echo  squall  lines  is  the  erhanced  mid-level  inflow  on  the  storms 
upshear  flank.  The  sequence  of  events  in  the  2  August  squall  line 
suggest  that  this,  feature  developed  in  response  to  a  strong  downdraft 
and  thus  the  severe  surface  winds  preceeded  the  development  cf  a  strong 
rear  inflow  jet  at  mid-levels.  The  question  remains  whether  this  is  a 
typical  or  atypical  sequence  of  events  and  as  such  will  be  the  subject 
of  fui  ther  inquiry.  Fina11y  the  dynamics  of  storm,  systems  such,  as  bow- 
echo  squall  lines  and  other  subclasses  of  convective  organization  that 
may  form  ir.  an  ensemble  of  convective  cells  (particularly  these 
conducive  to  severe  wind  events'  need  to  be  better  defir.ee  and  hopefully 
predicted.  7*.  is  hoped  that  these  rrief  reports  provide  some  guidance 
for  future  research  ar.d  understanding. 

5.0  5UMMAF  Y  CF  FA  Ki  CFYELCPMFNT 

PAMS  has  beer,  undergoing  a  major  overhaul  and  enhancement  to  allow 
several  new  features  to  be  added  to.  ’he  system  ar.d  tc  correct  for  seme 
old  deficiencies.  Sore  reaso-t  for  the  majer  changes  ir.  the  system  are: 
®  to  allow  fc-r  the  implementation  cf  twe- way  nesting  capabilities  in  the 
model 

t  to  add  the  comp  utational  schemes  cf  the  "Fieike "  model  to  the  HAMS 
system 

«  to  allow  for  faster  execution  time  or,  computers  which  require  longer 
vector  lengths 

®  to  conform  as  closely  to  the  FORTRAN  77  standard  as  possible  to 
facilitate  portability 

Tc  accomplish  these  goals,  sever;il  changes  needed  to  be  made: 

•  A  new  preprocessor  was  written.  Itself  in  the  FORTRAN  77  standard,  to 
correct  many  of  the  problems  with  the  old  PP1 .  Obsolete  features  of 


PP1  were  removed  (e.g.,  IF  Hock  processing  which  KFTFAN  77  car. 
handle).  The  new  preprocessor  also  is  much  faster  executing,  taking 
or.ly  about  1/6  as  long  as  PP1  . 

a  The  ole  FA  PS  was  run  with  a  ’’slat  "  I/O  scheme  where  all  the 
variables  cf  a  subset  cf  the  coiasir  was  kept  in  central  memory. 
Because  of  this,  most  cf  the  orc-C'  1  was  vectorio.ee  only  in  vertical 
loops .  These  short  loops  were  relatively  efficient  cr.  a  CFAY,  but 
were  much  tcc-  short  to  be  efficient  or.  a  CYfcEF  205 .  The  I/O 
structure  fer  the  r.ew  HAMS-  was  changed  tc  keep  the  whole  domain  ir. 
memory  but  only  a  subset  cf  the  model's  3-1  varieties.  This  scheme 
«.I  lews  full  3  T  '.reps  t '  be  uses  ir.  the  code  whf  r  e  a  CYPLP  oar 
vec*  cruse  the  .-cce  across  a.,  thr  le  ops.  This  also  eliminates:  the 
r.ect  ss:  ?  y  of  ’!>  macro  ca;  assist. ■:  cf  PPj,  making  devel  oprer.t  cf  tr. 
code  much  easier. 

0  The  new  IT  •??  rust  ur  c-  cf  th-  r.cw  '-AM-'  make  r  t  he  ir":  ere  r.taiic  r.  cf  tc 
the  two-  wav  nesting  a  1  g;  r  i  thus  the  tuslviicr  r.  f  the  "Fie  Ike  " 

schemes  much  caste:.  For  the  r.c  4  .  ;.g  it  :s  ue  c tc  inter  pt  lat 
tr.e  t<  ur.de  i  y  .v  ::  ti  tiers  of  the  r.ests  free  the  ec  a  rs.tr  rests. 

Tr.if  prccess  is  r  .h  easier  if  tr.'.  whole  ter- air.  cf  a  vtr.able  :s  :.r. 

r-f  a*  (.St  .  her  the  "F’.eihc"  s  or  ernes ,  cash  v r :  a  1 1  e  ; .  u;  : 
after  evert  terr.  c  a  1  r -ul  c  tec .  “r  tr  would  hast  beer,  economical  ly 

::f‘  sir  I-  ::  ‘it  slab  l.1  scheme.  Tr.  adci  tic-r,  there  art  :r_e 
:  ..mer:  sa!  :  •  :  iier  which  reqvirt.  the  entire  domain  in  •  c.re  <c  .g.  .  .  uh 

ss’.int  aeve-.  t  ior. ,  1  ■  •»  ' :  filter).  1  t-ngthy  sorat'.h  -trays  v:  IS  have 

see...  ncedeu  w  i  t  h  '  he  -1.  ;  I/O  r *  r.  -  . 

ft  With.  the.  ctr;:  n:  to  serf  err.  tc  the  l-'C  RTF  Ah  77  stardaro,  many  cf  the 

individual  p  by;  ir  al  cor:;  ;t«.iti  ur.ai  re  i  sues  such  as  the  it  it  r  c  p  hyri  >  s 


and  radiation  par ane t er iza tior.s  nave  been  rewritten  in  a  very  general 
form  so  that  they  could  be  transported  to  even  a  different  model  in 
the  future  if  desired.  A  summary  of  RAM?  is  given  in  Appendix  2. 

Test  « i m:.:] £ t i ons  have  recently  begun  with  the  new  svsten. 
Comparison  runs  with  the  old  "Pielke”  model  have  had  favorable 
results  with  only  a  few  minor  discr epanoies  remaining.  A  two-way 
nested  sea  breeze  simulation,  both  hydrostatic  and  no n- hydro s ta t i c- . 
has  beer  run  to  exercise  the  nesting  algorithms  and  again  reasonable 
results  were  produced.  It  is  anticipated  that  a  preliminary 
''production"  version  of  the  r.cv  RAMS  will  be  available  within  a  few 
months.  A  preprint  of  a  paper  presented  at  the  Er.viroscft  ’86 
conference  in  ten  Angeles  in  November  1986  which  gives  further  details 
of  the  new  system  and  it?  options,  is  reproduced  a.~  Appendix  3. 

The  hydrostatic  rested  re, ler?l  scale  version  cf  RAM?  is  being 
evaluated  ty  simulating  rerrs'-ale  convective  complexes  observed  during 
the  r '’0 F F  and  f FACF  field  programs.  The  simulations  are  being  rur. 
using  version  6  of  RAMS,  the  hydro-static  model  with  forward- backward 
‘  :  me  differ  er.;  ir.g.  Farr,  fa.'-  :  r  set.  up  :•  irr.llar!  y,  £  coarse  grid  is 
rur  with  about  r  lit'  kn  rercl  nt :  cr.  and  a  Kuo-type  convective 
parameter inatior  over  a  dc-m;  :  of  about  4500  kc  .  A  finer  grid  cf 
about  40  km  is  ther  one-way  nested  within  the  coarser  grid.  A 
Frits ch.-Chappell  type  convective  parameterization  is  used  or.  this 
higher  resolution  grid.  Sensitivity  experiments  are  rur.  or.  the  finer 
grid  testing  the-  sensitivity  to  various  physical  and  numerical 
processes.  From  the  results  so  far,  several  statements  can  be 
inferred  from  the  model  results. 

•  The  development  of  the  mere  seal  *•  circulation  of  the  convective 


system  is  dependent  cn  the  convective  heating  which,  in  turn,  is 
dependent  on  many  of  the  more  arbitrary  aspects  of  the  convective 
parameterizaticns. 

®  For  the  development  of  the  mesoscale  anvil  circulation,  the 
convective  heating  must  be  widespread  and  persistent,  not 
necessarily  very  strong.  The  heating  must  be  projected  or,  a  scale 
larger  than  the  Rossty  radius  of  deformation  and  over  a  longer 
timescale  than  that  of  an  individual  cumulus  cloud  to  form  a  quasi- 
balanced  circulation. 

®  For  the  convection  to  be  that  lor.g-iived  and  widespread,  the 
synoptic  and  mesc-a  scales  must,  provide  support  (conditional 
instability,  low-level  convergence,  upper- level  divergence,  etc.) 
for  the  continued  convection. 

«-  The  inclusion  of  a  resolvable  scale  ice  phase  ar.d  the  associated 
increase  ir.  the  amount  cf  latent  heat  releasee  helps  to  "spin  up" 
the  resolvable  scale  circulation  ar.d  increases  the  magnitude  of  the 
mesoscale  circulation. 

®  The  propagation  of  the  simulated  f-'CC  ic  dependent  or.  the  behavior  c 
the  low-level  jet,  even  ir.  this  tarociinic  environment.  As  the 
Right  progressed  ar.c.  the  jet  veered,  the  system  moved  toward  the 
east  (at  approximately  the  observed  speed)  along  with  the  maximum 
convergence  zone  at  the  rc-s e  of  the  simulated  jet. 

e  The  timing,  location,  and  strength  cf  the  modelled  convective 

systems  is  sensitive-  to  the  model  initial  conditions,  especially  to 
the-  initial  state  of  the  soil  moisture  (which  is  r.ct  routinely 
observed!.  Variations  in  the  soil  moisture  car.  dramatically  affect 
the  PbL  structure  and  therefore  affect  the  atmospheric  circulation 


in  the  vicinity  of  gradients  in  PBL  characteristics. 

Other  simulations,  testing  the  sensitivity  of  the  convective 
systems  to  radiation,  microphysical  parameters,  and  precipitation  wiJ 
be  run  ir.  the  near  future. 
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ABSTRACT 

In  this  paper  we  present  some  analytical  calculations  for  Level  2.5 
and  3  of  the  well-known  Mellor  and  Yamada  (M-Y)  hierarchy  of  second- 
order  turbulence  closure  models.  For  Level  2.5  we  extend  previous 
calculations  by  M-Y  and  show  that  the  requirement  of  positive¬ 
definiteness  of  the  u"2,  ^  variances  imposes  stringent  conditions  on  the 

possible  values  of  mean  flow  quantities  G.,  and  G„  (non-dimensional  mean 

r  H 

flow  shear  and  temperature  gradient,  respectively) . 

Two  additional  parameters  that  are  related  to  the  predictive  flow 

2  2  2 

2  2  -  -  ~ 

variables  G  =  “  and  a  =  (*p)  /[(— ')  +  (*p-)  ]  are  shown  to  be 

122  dz  dz  dz 

Q  Q 

important  to  describe  the  eddy-exchange  coefficients  and  fluxes  for  the 
Level  3.0  and  2.5  model.  Some  of  the  Level  3.0  functions  are  presented, 
and  the  possible  range  of  model  applicability  is  indicated. 

Finally  we  comment  on  the  problem  of  realizability  conditions  and 
applicability  of  the  Mellor  and  Yamada  models. 
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1.  INTRODUCTION 

In  a  sequence  of  papers,  Mellor  and  Yamada  (Mel lor , 1973 , 1 97  7 ; 

Mellor  and  Yamada, 1974, 1982;  Yamada , 1975 , 1983 ;  Yamada  and  Mellor, 1979) 
(hereafter  refered  to  as  M-Y)  developed  a  hierarchy  of  models  based  on 
the  second-order  closure  approach.  Different  degrees  of  simplification 
were  employed,  and  the  resulting  schemes  were  defined  as  levels . 
Initially,  four  levels  were  introduced  (Level  4  the  most  complex  and 
Level  1  the  simplest) ,  but  soon  it  became  apparent  that  this 
classification  was  too  narrow,  and  now  we  have  levels  1.0,  1.5,  2.0, 
2.5,  3.0,  3.5  and  4.  The  M-Y  hierarchy  of  models  gained  significant 
popularity  and  was  applied  by  many  to  various  atmospheric  flows  (see 
Mellor  and  Yamada  1982  review  paper  for  a  virtually  complete  list  of 
their  hierarchy  of  models  applications). 

It  was  recognized  by  M-Y  (1982)  and  others  (Hassid  and  Galperin, 
1983;  Helfand  and  Labraga,  1985)  that  for  some  combination  of  the  non- 
dimensional  mean  flow  parameters  (these  can  be  obtained  from  the 
prognostic  equations)  the  model  predicts  negative  variances  or  results 
in  singularities  in  solutions.  Such  non-physical  behavior  appears  to  be 
related  to  highly  transitory  phenomena  such  as  rapid  change  from  day  to 
night  conditions  or  to  the  rapid  entrainment  at  the  top  of  the  mixed 
layer . 

Mellor  and  Yamada  mention  that  negative  variances  appear 
occasionally  in  their  simulations.  In  the  1974  paper  they  state  that: 
"The  components  u^,  v^,  and  w^  should,  of  course,  be  positive  definite. 
This  turns  out  to  be  true  in  practice  everywhere  in  the  Level  2 
calculation  and  nearly  so  for  Levels  3  and  4  where,  however,  small 
negative  values  appeared  between  0700-0800  (after  a  discontinuity  in  the 


1 00 


tendency  of  wall  temperature)  at  a  couple  grid  points.” 

Recognition  of  such  non-physical  behavior  lead  M-Y  and  others 

(Hassid,  1983;  Helfand,  1985)  to  the  so-called  realizability  conditions 

for  possible  values  of  turbulent  fluxes  and/or  mean  flow  quantities. 

Simply  stated,  this  means  that  whenever  a  model  gives  non-physical 

behavior,  the  clipping  approximation  (setting  a  non-physical  value  to 

some  other,  usually  arbitrary,  value)  is  employed  or  change  of 

(predicted)  mean  flow  gradient  has  to  be  made.  Therefore,  it  seems  to 

be  of  importance  to  predict  the  possible  limits  of  M-Y  hierarchy  of 

models  in  advance.  This  provides  insight  on  the  range  of  model 

applicability  and  saves  one  frustration  during  the  code  developement . 

In  Chapter  2  we  discuss  briefly  M-Y  Level  2.5  and  3.0  models  and 

resulting  eddy  exchange  coefficients  and  turbulent  fluxes.  In  Chapter  3 

we  present  some  analytical  calculations  for  the  Level  2.5  and  3  models. 

For  Level  2.5  we  extend  previous  calculations  by  M-Y  (1982)  and  show 

that  the  requirement  of  positive-definiteness  of  the  u^  and  v^  variances 

imposes  stringent  conditions  on  the  possible  values  of  the  mean  flow 

quantities  G  and  G..  (non-dimensional  mean  flow  shear  and  temperature 
M  H 

gradient,  respectively). 

Two  additional  parameters  that  are  related  to  the  predictive  flow 

2  2  2 

2  2  -  -  -  z 

variables  G  =  ~ -  and  a  =  (fp)  /K^P)  +  (*p)  ]  are  shown  to  be 

T  2  2  dz  dz  dz 

Q  q 

important  to  describe  eddy-exchange  coefficients  and  fluxes  for  the 
Level  3.0  and  2.5  model.  Some  of  the  Level  3.0  functions  are  presented 
and  possible  range  of  model  applicability  is  indicated. 

Finally,  in  Chapter  4,  we  comment  on  the  problem  of  realizability 
conditions  and  applicability  of  the  Mel  lor  and  Yamada  models. 
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2.  MODEL  EQUATIONS 

M-Y  use  departure  from  isotropy  as  a  small  parameter  in  their 
analysis  of  the  Reynolds-averaged  equations.  They  show  that  in  the 
atmospheric  boundary  layer  (ABL)  one  can  neglect  diffusion,  time  and 
adveetion  terms.  Level  3  predicts  only  variances.  In  Level  2.5  only 
the  equation  for  the  turbulent  kinetic  energy  (TKE)  is  retained. 

Because  of  simplifications  applied  to  the  theory,  the  resulting  set 
of  Reynolds  equations  is  linear  and  can  be  solved  analytically.  The 
algebra  used  in  this  study  is  extensive  (see  Yamada ,  1978;  Flatau,  1985 
for  more  details).  For  many  algebraic  manipulations  in  this  paper  we 
used  the  symbolic  algebra  program  REDUCE  3  (Hearn, 1984) ,  We  will 
present  results  for  the  calculations  of  the  dry  case  (i.e.,  without 
cloud  microphyslcs-turbulence  covariances).  The  extension  to  a  moist 
c«se  is  relatively  simple. 

To  ma^e  the  final  formulas  more  compact,  let  us  introduce  several 

non-dimensional  quantities.  The  non-dimensional  eddy  exchange 

coefficients  S..  and  S„  are  defined  by  the  following  equations 
M  H 

Km  =  lqSM  (1) 


and 


Kh  -  lqSH 


(2) 


where  K..  and  K.,  are  momentum  and  heat  eddy  exchange  coefficients,  1  is 
M  H 

turbulent  length  scale,  and  q  is  defined  by  q2  =  u^  +  v^  +  w^  .  The  uS 
^  and  ^  are  velocity  variances. 


We  also  define  the  non-dimensional  parameters  Gu,  G„,  G„,  and  a, 

M  H  I 


which  are  related  to  the  mean  flow  variables 


cM  ■  ^  (f  )2  * 

H  2  dz  dz 

q 


r  1  .  a© 

GH  =  “~2  **& 

q 


G  ,  li  U&l-al 

°T  2  2  S 

q  q 


(f>  /[(f)  ♦  (f )  ] 

dz  dz  dz 


where  0  =  1/0  ,  ©  is  the  basic  state  potential  temperature,  6Z  is  the 

oo 

du  av 

potential  temperature  variance,  —  ,—  are  mean  flow  velocity  gradients, 

a© 

and  7-  is  the  potential  temparature  gradient.  It  can  be  seen  that 
0  z 

0  <  a  <  1 .  Notice  that  G  ,  G  ,  G  ,  and  o  are  related  to  the  prognostic 

H  I 

equations  in  the  model.  In  other  words,  if  we  express  all  variables  in 
terms  of  G^,  G^ ,  G^,  and  a,  then  the  problem  is  solved.  After  some 
algebra  we  can  show  that  S  and  S„  for  the  Level  3.0  are  given  by 

sn 

s  - 

SK  -  cd' 


S" =  Va 


(7a, b) 


where 

SM  =  9GHA2A1(A2  +  4A1C1)  +  a1(27AjGt  +  36A2AjGt  -  3C  ♦  1) 


SH  •  -^*2*1  *  180hV2‘!C1  4  W^Wt  *  11  *  18V2‘!GT  *  3A2°T 
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r, 

L 

L 


Sd  ‘  108C5*2*1  *  5«HGH*2*»  -  21GB#2»1  *  ‘Vl  *  1 

These  formulas  are  the  same  as  those  in  Mellor  and  Yamada's  (1974)  paper 
(eqs  55,56),  although  they  are  written  here  in  the  non-dimensional  form. 

For  the  Level  2.5  we  get  the  non-dimensional  eddy  exchange 
coefficients  in  the  form 


S 


M 


S 


H 


(8a,b) 


where 

Sd  '  *  120hVi  -  11  "V? 

-  9GHA2V~GMGH(MA2A1  ♦  72A2aJ)  +  G^OE,,  +  12A1)  "  1 

■  -<180^  -  90^^  .  1)A2 

SK  ’  -Vl<9CIB2*2  *  36C1*2,1  •  3E2*2  *  9A2  ‘  C,  *  11 

For  both  levels  the  final  set  of  equations  is 

^/q2  =  1/3  ♦  Al6SMGM  -  2AlSHGH, 

^/q2  =  1/3  +  -  2AlSHGH.  and 

w^/q2  =  1/3  ♦  2A:SmGm  ♦  4AlSHGH, 


where  6  =  2(3a  -  1)  and  8  =  2(2  -  3a)  . 


—  ..  . 2  3U8V 

uv  =  6A, 1  Su  r  V" 
1  M  dzdz 


uw  =  - 1 q  S 


3U 
M  dz 
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vw  =  -lqS 


3V 
M  3z 


“e  -  3V,SH  *  V  51 

55  -  3A21(Ss  *  SM»  ff 

—  3© 

wG  = 


dz 


av 

0Z 


and  for  Level  2.5  only 


51  -  v3„f! 


( 1  Oa .  b ,  c  ' 


( 1 1  a ,  b ,  c ) 


(12.' 


The  ratio  of  TFCE  production  to  dissipation  can  be  expressed  as 

<Ps  *  V/e  ■  E1,SMGM  *  W  "3I 

3.  ANALYTICAL  RESULTS  FOR  THE  MEXLCR  AND  YAHADA  HIERARCHY  OF  MODELS 

In  oiiis  chapter  we  discuss  some  analytical  results  and  realizabilty 
conditions  for  the  Mellor  and  Yamada  (1982)  hierarchy  of  models  -  Levels 
2.3  and  3.0.  The  purpose  of  such  a  study  is  to  decide  what  the  possible 
limitations  of  the  scheme  are.  Using  the  typical  values  of  parameters 
G^,  G^,  G^,,  Rj,,  and  a  (see  Chapter  2).  we  are  able  to  plot  fields  of 
turbulent  variances  and  co-variances .  Obvious  criteria  such  as 
positive-definiteness  of  variances  are  used  to  limit  the  possible  range 
of  mean  flow  values. 

3.1  Level  2.5 

The  Level  2.5  model  depends  on  three  predictive  parameters.  These 
are  G^  (non-dimensional  temperature  gradient) ,  G^  (non-dimensional 
quantity  related  to  the  mean  flow  shear),  and  a  (ratio  of  shear  squares 
in  the  x  and  y  directions).  See  Chapter  2  for  their  definitions.  The 
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contours  nf  turbulent  quantities  such  as  non-dimensional  momentum  ar.d 

heat  eddy  coefficients,  turbulent  variances  u^  ,  and  ratio  of 

tr.e  turbulent  production  to  the  turbulent  dissipation  are  plotted  here 

with  G0  and  G..  as  independent  variables.  The  a-dependence  of  these 
H  M 

turbulent  quantities  will  also  be  discussed. 

Let  us  start  with  the  Figure  2b  where  (P  +  P.)/e  is  plotted.  The 

s  b 

G..  values  are  always  positive  and  are  shown  there  for  (0.,3.)  region. 

The  G  values  can  be  negative  ( for  stable  case)  and  positive  (for 
H 

unstable  uase).  The  (-0.4, 0.1)  range  of  G„  values  is  used  here.  The 

ri 

isoline  (P  +  P  )/e  =  1  corresponds  to  the  local  equilibrium  assumption 

5  D 

(i.e.,  it  satisfies  the  Level  2.0  model  conditions).  On  the  other  hand, 

for  arbitrary  G,,  and  G.,,  the  (P  +  F.  )/e  4  1.  This  is  clearly  seer,  in 
h  M  s  b 

Figure  2b.  In  Figure  2b  isolir.es  labeled  INF-1  represent  values  of  G(, 

ar.d  G..  for  which  ‘  P  ♦  P  )/e  singular.  Curves  labeled  PI  and  P2 
M  s  b 

give  G  and  G  values  where  (F  *  FK)/e  8°es  to  zero.  Contours  are 

n  rl  S  u 

plotted  for  (P  -  Pu)/e  -  ( 0 . 5 . 1 . 0 , 1 . 5 . 2 . 0 , 2 . 5 ) . 
s  b 

The  PI  isoline  originates  at  Gw  =  0  and  G,.  =  0.  It  delimits  the 

K  h 

regicn  of  positive  (upper  part)  and  negative  (lover  part)  values  cf 

(P  +  P,)/e.  For  the  stable  case  (G„  <  0)  and  small  Gw,  the  buoyancy 
s  b  H  M 

term  P^  becomes  negative  and  is  larger  than  shear  generation  term  P^ . 

This  leads  to  negative  (P  +  F  )/e.  On  the  other  hand,  for  the  unstable 

s  b 

case,  values  of  (P^  +  are  always  positive. 

Figures  1  and  2a  present  values  of  eddy  exchange  coefficients  S 

n 

and  S^j.  The  INF-1  curve  is  the  same  as  before  and  represents  G^.  and  G^ 
values  for  which  S„  and  Sw  tends  to  infinity.  These  large  values  of 

n  M 

exchange  coefficients  are  not  physically  possible,  and  the  region  close 


to  the  INF  curve  has  to  be  somehow  excluded.  The  values  to  the  left  of 


INF  curve  are  positive,  which  indicates  positive  exchange  coefficients. 

On  the  S..  diagram  the  SM  curve  is  for  S  =  0.  In  the  region  between  INF 
M  M 

and  SM  curves  in  Figure  1,  the  S  values  are  negative:  however,  from 

.  . 

other  considerations  (positive  definiteness  of  velocity  variances),  we 

will  have  to  exclude  this  region.  In  other  words.  Level  2.5  is  forced 

to  predict  positive  (down-gradient)  values  of  exchange  coefficients. 

Positive  definiteness  of  velocity  variance  also  leads  to  some  limits  or. 

the  possible  values  of  G..  and  G... 

H  M 

Compare  now  Figures  1  (S  )  and  2b  (production/dissipation) .  Notice 

M 

that  the  gradient  Richardson  number  is  giver,  by  Ri  =  -G ,/G^..  One  car. 

H  M 

then  easily  obtain  Sw  as  a  function  of  Ri  for  constant  values  of 

r. 

(P  +  P. )/e.  Figure  6  combines  plots  of  the  curves  discussed  above 
s  b 

(i.c.,  S  ,  S  ,  and  (P  +  P,  )/t  )  and  can  be  used  for  tracing  S  ,  Su 
M  H  s  b  H  M 

behavior  for  constant  ( P^  +  P^J/e  values.  In  Figure  6  heavy  solid  lines 

are  S  values,  S  values  are  represented  by  broken  lines,  and  thir.  solid 
M  H 

lines  are  (P  +  P  )/e  values.  For  G  approaching  0  with  the  constraint 

S  D  M 

(P  -1-  P,_ )  /  e  =  const.  (i.e.,  traveling  down  the  (P  *  F  )/e  isoiir.e),  we 
s  b  s  b 

get  a  finite  value  of  S^.  This  limited  value  corresponds  to  infinite 

(and  negative)  Richardson  number.  Notice  that  (P  +  P  )/e  >  0  isolir.es 

s  b 

cross  the  G  axis  for  positive,  non -zero  values  of  G  .  The  above 
r.  H 

discussion  applies  to  the  other  turbulent  covariances.  They  all  can  be 
expressed  as  a  function  of  only  one  mean  flow  quantity  (Richardson 
number)  if  a  certain  constraint  on  the  ratio  of  the  turbulent 
production  to  the  turbulent  dissipation  is  imposed. 

Figures  3-5a,b  show  turbulent  variances  of  velocity  components 
wz/q  ,  uz/q  ,  vz/q  .  From  Eqs.  9a, b  we  see  that  u-  and  v-variances 
depend  on  a.  It  can  be  shown  that  from  the  point  of  view  of 
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realizability  conditions  the  a  =  0  and  a  =  1  give  the  most  stringent 

limitations  on  the  possible  G  values.  Therefore,  all  the  results  here 

H 

are  presented  for  a  =  1;  i.e.,  for  ~  =  0.  The  w^/q^  is  not  dependent 

— T  2 

on  a.  Examining  Figures  3-5  we  see  that  v^/q  is  positive  onlj  for  Gu  < 

H 

0.0288,  and  this  limit  seems  to  guarantee  positive  definiteness  of  all 
velocity  variances. 

Finally  we  present  velocity  variances  and  (F  +  P.)/e  as  a  function 

s  b 

of  S  and  S  (Figures  7-8a,b).  Again  a  =  1.  Only  positive  isolines  and 

r.  ri 

those  that  are  smaller  than  one  are  plotted  for  the  velocity  variances. 

The  isolines  are  labelled  A.B.C.D,  and  E  and  correspond  to  (0.001, 

0.2,  0.5,  0.7,  1.0)  values.  The  isotropic  case  (all  variances  equal  to 

1/3)  gives  approximately  the  same  values  of  S  and  S^.  The  realizable 

solution  follows  this  pattern  for  other  velocity  variances;  i.e.,  the  S^_, 

and  S..  values  are  of  the  same  order, 
n 

3.2  Level  3.0 

This  level  of  simplification  was  considerably  less  investigated  and 
used  by  researchers.  We  performed  analysis  of  this  model  and  some  of 
the  results  are  presented  in  Figures  9-14.  Here  the  additional 
parameter  G„  comes  to  play,  which  makes  analysis  more  involved.  We  plot 
G  -G„  dependence  of  turbulent  quantities  for  2  different  values  of  G„  = 

r.  n  i 

(0.1,  3.0).  The  range  of  Gu  and  G„  is  from  -6  to  6  on  all  plots. 

n  n 

Values  of  G^  have  to  be  positive,  but  negative  ones  are  retained  here 

for  illustrative  purposes.  As  before,  negative  values  are  plotted  as 

2 

dashed  lines.  Contours  for  the  normalized  (by  o  )  velocity  variances 
are  chosen  at  0.05,  0.2,  0.5  and  1.0  (0.12  instead  of  0.2  for  the 
vertical  velocity  variance).  By  definition,  all  velocity  variances  have 
to  be  in  the  (0,1)  range.  One  can  see  (e.g.,  in  Figure  9)  that  values 


giving  realizable  solutions  for  the  Level  2.5  may  now  correspond 

M  H 

to  negative  variances.  The  solutions  are  highly  variable  with  G  .  The 

realizable  values  are,  in  general,  located  close  to  the 

(P  +  P. )/e  =  1  isoline, 
s  b 

We  diagnostically  calculated  (Flatau.  1985)  Level  3.0  functions 
from  the  numerical  results  of  our  runs  of  the  Level  2.5  model.  The 
results  indicated  that  the  Level  3.0  produces  very  stringent  conditions 
on  the  G,,  and  Gu  values.  The  picture  is  still  not  clear;  however,  we 

rJ  H 

hope  to  obtain  some  diagnostic  results  from  the  working  Level  3.0  model. 
4.  CONCLUSIONS 

As  was  mentioned,  it  seems  that  the  Level  2.5  and  3.0  model  does  not 
properly  handle  transitory  phenomena  but  has  the  ability  to  adjust  the 
mean  flow  and  turbulence  to  consistently  physical  values.  As  long  as 
variances  are  not  used  explicitly  (for  example  to  diagnose  triple 
correlation  terms)  but  are  just  diagnosed,  the  fact  that  they  are 
negative  is  only  annoying.  The  slowly  varying  phenomena,  e.g.,  dry 
boundary  layer  on  a  sunny  day  -  as  in  the  Wangara  case  -  are  predicted 
well  by  the  lower-level,  second-order  closure  schemes.  The  apparent 
failure  of  Level  3.0  and  2.5  to  predict  highly  variable  flows  can 
probably  be  traced  to  the  basic  assumption  of  Mellor  and  Yamada  - 
namely,  scaling  out  the  advection  ar.d  diffusion  terms  in  almost  all 
equations.  At  the  top  of  the  stratocumulus  layer  where  intensive, 
local  (radiative)  cooling  is  distributed  mainly  through  turbulent 
diffusion  and  where  gradients  of  stratification  are  strong,  the 
assumption  of  the  small  role  of  diffusional  terms  can  be  a  weak  one  (see 
e.g.  Chen,  1984  simulations  of  Sc  deck).  In  such  a  case,  when  greater 
definition  of  the  turbulence  is  important,  our  results  indicate  that  the 
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higher-level  second-order  turbulence  closure  scheme  should  be  used.  It 
is  true,  though,  that  the  model  is  driven  towards  local  equilibrium. 
Imagine  vigorous  TKE  advected  to  an  otherwise  undisturbed  region. 
Parameterization  of  dissipation  provides  negative  feedback  in  the  TKE 

equation 

,  2  3 

is_  =  _s_ 

at  ’ 

and  the  model  will  adjust  to  local  equilibrium.  This  fact  is  probably 
the  only  Justification  of  usage  of  the  Level  2  (which  is  based  on  the 
local  equilibrium  assumption)  model  together  with  a  prognostic  equation 
for  the  TKE  (this  scheme  was  used  by  Yamada  (1983)).  Such  a  scheme,  by 
definition,  cannot  take  into  account  situations  which  are  strongly  out 
of  equilibrium.  Similiary,  we  have  shown  in  Chapter  3  that  the  Level 
2.5  and  3.0  give  realizable  solutions  only  for  the  situations  close  to 
the  local  equilibrium. 

Probably  the  simplest  and  the  strongest  recommendation  that  can  be 
inferred  from  this  study  is  to  apply  the  lower-level  second-order 
schemes  in  accordance  with  their  assumptions.  The  Mellor  and  Yamada 
hierarchy  is  based  on  small  parameter  scale  analysis  (departure  from 
isotropy),  and  if  this  assumption  is  not  satisfied,  the  scheme  will  not 
work  properly.  Although  the  realizability  conditions  should  be  checked 
diagnostically,  their  use  should  be  limited  as  much  as  possible.  If  the 
model  is  consistently  not  realizeable,  the  only  reasonable  solution  is 
to  implement  the  higher-level  scheme.  If  the  non-real izable  solutions 
occur  sporadically  (model  spin-up  time,  transition  from  day  to  night, 
top  of  the  PBL) ,  and  if  the  main  emphasis  of  the  model  is  not  on  a  very 
detailed  physical  depiction  of  these  regimes,  it  can  be  used  with  some 


confidence . 
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Figure  1.  Level  2.5  functions.  Constants  of  the  model  are  indicated  at 
the  top  of  each  plot.  They  are  in  the  Aj ,  A2 ,  Bj ,  B2,  C1(  C2, 
C3  order.  The  INF-1  curve  Indicates  values  where  functions 
tend  to  infinity.  Here  momentum  eddy  exchange  coefficient  Sp 
is  plotted. 

Figure  2.  (a)  Same  as  in  Figure  1,  but  for  S^;  (b)  same  as  in  Figure  1, 

but  for  (P  +  P.  )/e  ,  PI  curve  indicate  points  where  function 
s  b 

goes  to  0. 

— T  2 

Figure  3.  (a)  Same  as  in  Figure  1,  but  for  wz/q  ,  VEL-1  is  where 

function  goes  to  0,  (b)  same  as  (a),  different  GH,  G^. 

Figure  4.  (a)  Same  as  in  Figure  1,  but  for  u^/q2,  VEL-1  is  where 

function  goes  to  0,  (b)  same  as  (a),  different  GH,  . 

Figure  5.  (a)  Same  as  in  Figure  1,  but  for  v^/q2,  VEL-1  is  where 

function  goes  to  0,  (b)  same  as  (a),  different  G^,  GM. 

Figure  6.  Composite  of  plots  presented  in  Figures  1  and  2.  The  heavy 
solid  lines  indicate  S^  functions,  broken  lines  are  isolines 
of  Sfj,  and  thin  solid  lines  are  isolines  of  the  ratio  of 
turbulent  production  and  dissipation.  The  INF-1  curve  gives 
singular  points  for  all  functions.  PI  and  P2  curves  are  for 


(P  +  P  ) / e  =  0. 
s  b 


Figure  7.  Level  2.5  plots  as  a  function  of  and  S^.  (a)  For  u^/q2, 

(b)  for  v^/q2. 

Figure  8.  Level  2.5  plots  as  a  function  of  and  (a)  w^/q2.  (b) 

(P  +  P  )/e . 
s  b 

Figure  9.  Level  3.0  functions.  Constants  of  the  model  are  indicated  at 
the  top  of  each  plot.  They  are  in  the  Aj ,  A2,  Bj ,  B2 ,  Cj,  Gy. 
The  INF-1  curve  indicates  values  where  functions  tend  to 
infinity.  Here  Gj  =  0.3.  Numbers  at  the  bottom  of  each  plot 
indicate  values  for  Gy  =  0  where  the  denominator  (first 
two  numbers)  and  numerator  (next  number  or  numbers)  of  the 
function  tend  to  zero. 

Figure  10.  (a)  Same  as  in  Figure  9,  but  for  (P  +  Pw)/e,  the  PI  curve 

s  b 

indicates  values  of  G^  and  G^  where  (Ps  +  P^J/e  goes  to 
zero,  (b)  same  as  (a)  but,  for  u^/q^. 

Figure  11.  (a)  Same  as  in  Figure  9,  but  for  v^/q2,  (b)  same  as  ir. 

Figure  9,  but  for  w^/q2. 

Figure  12.  Same  as  in  Figure  9,  but  for  Gj  =  3.0. 

Figure  13.  Same  as  in  Figure  10,  but  for  Gj  =  3.0. 
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Figure  14.  Same  as  in  Figure  11,  but  for  Gj  =  3.0. 
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INTRODUCTION 

The  numerical  atmospheric  models  developed  independently  under 
the  direction  of  William  R.  Cotton  and  Roger  A.  Pielke  have 
recently  been  combined  into  the  CSU  Regional  Atmospheric 
Modelling  System  (RAMS).  RAMS  is  a  general  and  flexible 
modelling  system  rather  than  a  single  purpose  model.  For 
example,  current  research  using  RAMS  includes  atmospheric 
scales  ranging  from  large  eddy  simulations  (Ax  s  100  m)  to 
mesoscale  simulations  of  convective  systems  (Ax  *  100  km) . 

This  paper  will  discuss  the  options  available  in  RAMS,  the 
engineering  aspects  of  the  system  and  how  the  flexibility  is 
attained. 

RAMS  OPTIONS 

RAMS  is  a  merging  of  basically  three  models  that  were  designed 
to  simulate  different  atmospheric  circulations.  These  were  a 
non-hydrostatlc  cloud  model  (Tripoli  and  Cotton1)  and  two 
hydrostatic  mesoscale  models  (Tremback,  et  al.^  and  Mahrer  and 
Pielke^).  Because  of  this,  the  modelling  system  contains  many 
options  for  various  physical  and  numerical  processes.  These 
options  are  listed  below. 

1)  Basic  equations: 

Option  1  -  non-hydrostatlc  time  split  compressible, 
(Tripoli  and  Cotton1). 

Option  2  -  hydrostatic  time  split  compressible  or 
incompressible,  (Tremback  et  al.2) . 

Option  3  -  quasi-hydrostatic  incompressible, 

(Song  et  al.4) . 

Option  4  -  hydrostatic  incompressible, 

(Mahrer  and  Pielke^). 
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2)  Dimensionality: 

1,  2,  or  3  spatial  dimensions. 

3)  Vertical  coordinate: 

Option  1  -  Cartesian 

Option  2  -  sigma  -  z 

4)  Horizontal  coordinate: 

Option  1  -  Cartesian 

Option  2  -  latitude/longitude 

5)  Grid  stagger: 

"Arakawa  -  C  grid" 

6)  Finite  differencing: 

Option  1  -  leapfrog  on  long  timestep,  forward-backward 
on  small  timestep.  2nd  or  4th  order  flux 
conservative  advection. 

Option  2  -  forward-bu'kward  time  split,  6th  order  flux 
conservative  (Tremback  et  al.5)  or  cubic 
spline  advection. 

7)  Turbulence  closure: 

Option  1  -  Smagorinsky  -  type  eddy  viscosity  with 
Richardson  number  dependence. 

Option  2  -  eddy  viscosity  as  a  function  of  a  prognostic 
turbulent  kinetic  energy. 

Option  3  -  eddy  viscosity  from  O'Brien  profile  function. 

8)  Condensation: 

Option  1  -  grid  points  fully  saturated  or  unsaturated. 

Option  2  -  Kuo  -  type  convective  parameterization. 

Option  3  -  Fritsch/Chappell  -  type  convective 
parameterization. 

Option  4  -  no  condensation. 

9)  Cloud  microphysics: 

Option  1  -  warm  rain  processes,  (Tripoli  and  Cotton6). 

Option  2  -  option  1  plus  pristine  crystals  and  graupel, 
(Cotton  et  al.O . 

Option  3  -  option  2  plus  aggregates. 

Option  4  -  no  microphysics. 
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10)  Radiation: 


a.  Shortwave: 

Option  1  -  parameterization  described  in  Chen  and 
Cotton8. 

Option  2  -  parameterization  described  in  Mahrer  and 
Pielke8  . 

b.  Longwave: 

Option  1  -  parameterization  described  in  Chen  and 
Cotton8 . 

Option  2  -  parameterization  described  in  Mahrer  and 
Pielke8  . 

11)  Lower  boundary 

Option  1  -  specified  surface  temperature  and  moisture 
values  or  fluxes  with  surface  layer 
similarity  parameterization. 

Option  2  -  surface  values  obtained  from  surface  energy 
balance  using  a  prognostic  soil  model, 
(Mahrer  and  Pielke8). 

Option  3  -  modified  form  of  option  2  with  prognostic 

surface  equations,  (Tremback  and  Kessler9). 

Option  4  -  option  2  or  3  with  vegetation 

parameterization,  (McCumber  and  Pielke10) . 

12)  Upper  boundary  conditions: 

Option  1  -  rigid  wall. 

Option  2  -  gravity  wave  radiative  condition. 

Option  3  -  Raylegh  function  layer. 

Option  4  -  prognostic  surface  pressure  (hydrostatic 
options) . 

13)  Lateral  boundary  conditions: 

Option  1  -  Klemp/Wilhelmson  radiative  boundary. 

Option  2  -  Orlanski  radiative  boundary. 

Option  3  -  Klemp/Lilly  radiative  boundary. 

Option  4  -  option  1,2,  or  3  with  a  "mesoscale 

compensation  region”  (Tripoli  and  Cotton1). 

Option  5  -  Perkey  and  Kreltzberg11  sponge  boundary 
condition  when  simulation  is  run  with 
non-horizontally  homogeneous  initial 
conditions. 

14)  Initialization: 

Option  1  -  horizontally  homogeneous. 

Option  2  -  option  1  plus  variations  to  force  cloud 
initiation. 

Option  3  -  NMC  data  and/or  soundings  objectively 
analyzed  on  isentropic  surfaces  and 


interpolated  to  the  model  grid. 

Option  4  -  NMC  data  interpolated  to  the  model  grid. 


As  one  can  see,  RAMS  is  a  versatile  modelling  system.  It  has 
been  applied  to  a  variety  of  weather  phenomena  including: 

a)  towering  cumuli  and  their  modification 

b)  tropical  and  mid-latitude  cumulonimbi 

c)  mountain  slope  and  valley  circulations 

d)  orographic  cloud  formation 

e)  marine  statocumulus 

f)  mountain  wave  flow 

g)  large  eddy  simulations  of  power  plant  plume  dispersal 

h)  mountain  generated  convective  systems 

i)  mid-latitude  convective  systems 

J)  mesoscale  transport  and  dispersion 

k)  low-level  Jet  simulations 

l)  atmosphere  -  soil/ vegetation  interactions 

m)  sea  breeze  simulations 

ENGINEERING  ASPECTS 

Because  of  the  large  number  of  options  in  RAMS,  the  structuring 
of  the  code  needs  to  be  carefully  considered.  This  section 
will  discuss  various  aspects  of  the  code  structure  of  the 
system. 

Pre-processor  The  code  of  RAMS  is  written  in  as  close  to  the 
FORTRAN  77  standard  as  possible.  However,  with  a  program  as 
large  as  this,  the  FORTRAN  standard  is  lacking  in  several 
features  such  as  global  PARAMETER  and  COMMON  statements  and 
conditional  compilation.  To  remedy  these  lnsuf ficiences,  the 
RAMS  code  takes  advantage  of  a  pre-processor  written  as  part  of 
the  RAMS  package.  This  pre-processor  then  produces  standard 
code  which  can  be  then  compiled  on  a  FORTRAN  77  compatible 
machine.  The  pre-processor  itself  is  written  in  the  77 
standard  so  that  the  package  as  a  whole  is  highly  portable.  It 
takes  full  advantage  of  the  character  features  of  FORTRAN  and 
has  executed  successfully  on  a  number  of  machines  including  a 
VAX.  CRAY-1,  CRAY-XMP,  and  CYBER  205  without  modification. 

Some  of  the  features  of  the  pre-processor  are  described  below: 

1)  By  including  a  character  in  the  first  column  of  a 

line  of  code,  that  line  can  be  "activated"  or 
"eliminated”  from  the  compile  file.  This  allows 
for  conditional  compilation  of  single  lines  or 
entire  sections  of  code. 

2)  A  pre-processor  variable  can  be  set  to  a  value.  This 

variable  can  then  be  used  in  other  expressions 
including  a  pre-processor  IF  o^  block  IF  to 
conditionally  set  other  pre-processor  variables. 
These  variables  also  can  be  converted  to  FORTRAN 
PARAMETER  statements  which  can  be  inserted  anywhere 
in  the  rest  of  the  code. 


3)  A  group  of  statements  can  be  delineated  as  a 

"global"  which  then  can  be  inserted  anywhere  in  the 
code.  This  is  very  useful  for  groups  of  COMMON  and 
PARAMETER  statements. 

4)  DO  loops  can  be  constructed  in  a  DO/ENDDO  syntax. 

eliminating  the  need  for  statement  labels  on  the  DO 
loops. 

I/O  structure  For  those  machines  with  limited  central  memory 
and  a  "non-virtual "  operating  system  or  for  efficiency  on 
virtual  systems.  RAMS  is  constructed  with  a  disk  I/O  scheme. 
When  the  scheme  is  operating,  a  subset  of  the  model's  three- 
dimensional  variables  will  reside  in  central  memory  at  any  one 
time.  Computations  then  can  be  performed  with  this  subset. 

When  these  computations  are  finished,  a  new  subset  of  three- 
dimensional  variables  are  requested  and  computations  performed 
with  these.  The  RAMS  structure,  thus,  is  dependent  on  this  I/O 
scheme  and  consists  of  a  series  of  calls  to  the  I/O  scheme  then 
the  routines  which  do  the  calculations. 

Modularity  For  flexibility,  RAMS  is  written  as  modular  as 
possible.  Each  individual  physical  parameterization  or 
numerical  process  is  put  in  a  separate  subroutine  so  that  the 
routines  can  easily  be  replaced  for  different  options  or  with 
new  developments. 

Computational  routines  The  routines  that  do  the  actual 
computations  for  the  model  are  written  so  that  the  implementor 
of  a  new  or  replacement  routine  does  not  need  to  be  concerned 
with  most  of  the  details  of  the  rest  of  the  model  computations. 
All  three-dimensional  variables  are  "passed"  to  the 
subroutines  through  the  call  statement  with  other  variables 
passed  through  COMMON.  The  implementor  then  has  the 
flexibility  to  structure  his  routine  in  whatever  manner  he 
wishes  to  produce  the  desired  result.  This  concept  will  also 
make  the  implementation  of  routines  from  other  models  and 
programs  easier  with  less  modification  required. 

CONCLUSION 

RAMS  has  been  developed  at  Colorado  State  University  to  provide 
a  general  and  flexible  framework  for  many  scales  of  atmosheric 
flow  simulation.  A  pre-processor  is  used  to  extend  the 
capabilities  of  standard  FORTRAN  while  the  structure  of  the 
model  allows  the  implementor  much  freedom  in  the  developing  of 
individual  portions  of  the  model. 

RAMS  is  a  research  tool  and,  as  such,  is  not  a. static 
program.  Some  of  the  developments  planned  for  the  system  in 
the  near  future  are  two-way  interaction  nesting  capabilities 
and  sixth  order  leapfrog  advection. 
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